NASA/CR— 2000-210062 


ARL-CR-451 



Simulating Fatigue Crack Growth 
in Spiral Bevel Gears 


Lisa E. Spievak, Paul A. Wawrzynek, and Anthony R. Ingraffea 
Cornell University, Ithaca, New York 


May 2000 



The NASA STI Program Office ... in Profile 


Since its founding, NASA has been dedicated to 
the advancement of aeronautics and space 
science. The NASA Scientific and Technical 
Information (STI) Program Office plays a key part 
in helping NASA maintain this important role. 

The NASA STI Program Office is operated by 
Langley Research Center, the Lead Center for 
NASA's scientific and technical information. The 
NASA STI Program Office provides access to the 
NASA STI Database, the largest collection of 
aeronautical and space science STI in the world. 
The Program Office is also NASA's institutional 
mechanism for disseminating the results of its 
research and development activities. These results 
are published by NASA in the NASA STI Report 
Series, which includes the following report types: 

• TECHNICAL PUBLICATION. Reports of 
completed research or a major significant 
phase of research that present the results of 
NASA programs and include extensive data 
or theoretical analysis. Includes compilations 
of significant scientific and technical data and 
information deemed to be of continuing 
reference value. NASA's counterpart of peer- 
reviewed formal professional papers but 

has less stringent limitations on manuscript 
length and extent of graphic presentations. 

• TECHNICAL MEMORANDUM. Scientific 
and technical findings that are preliminary or 
of specialized interest, e.g., quick release 
reports, working papers, and bibliographies 
that contain minimal annotation. Does not 
contain extensive analysis. 

• CONTRACTOR REPORT. Scientific and 
technical findings by NASA-sponsored 
contractors and grantees. 


• CONFERENCE PUBLICATION. Collected 
papers from scientific and technical 
conferences, symposia, seminars, or other 
meetings sponsored or cosponsored by 
NASA.. 

• SPECIAL PUBLICATION. Scientific, 
technical, or historical information from 
NASA programs, projects, and missions, 
often concerned with subjects having 
substantial public interest. 

• TECHNICAL TRANSLATION. English- 
language translations of foreign scientific 
and technical material pertinent to NASA's 
mission. 

Specialized services that complement the STI 
Program Office's diverse offerings include 
creating custom thesauri, building customized 
data bases, organizing and publishing research 
results . . . even providing videos. 

For more information about the NASA STI 
Program Office, see the following: 

• Access the NASA STI Program Home Page 
at http:/lwiviv,$tLnasa,gov 

• E-mail your question via the Internet to 
help@sti.nasa.gov 

• Fax your question to the NASA Access 
Help Desk at (301) 621-0134 

• Telephone the NASA Access Help Desk at 
(301)621-0390 

• Write to: 

NASA Access Help Desk 
NASA Center for AeroSpace Information 
7121 Standard Drive 
Hanover, MD 21076 


NASA/CR— 2000-210062 


ARL-CR-451 



Simulating Fatigue Crack Growth 
in Spiral Bevel Gears 

Lisa E. Spievak, Paul A. Wawrzynek, and Anthony R. Ingraffea 
Cornell University, Ithaca, New York 


Prepared under Grant NAG3-1993 


National Aeronautics and 
Space Administration 


Glenn Research Center 


May 2000 



Acknowledgments 


The research contained in this thesis was conducted under grant NAG3-1993 between Cornell University and 
NASA Glenn Research Center. I wish to thank Dr. David Lewicki and Dr. Robert Handschuh of the U.S. 
Army Research Laboratory at NASA Glenn Research Center. Much of this thesis' work is a direct 
result of their advice and expertise. Lehigh University professor Dr. Eric Kaufmann's time and 
technical knowledge were instrumental with the scanning electron microscope observations 
contained in this thesis. In addition. Dr. Richard N. White at Cornell University 
volunteered his time and skills to photograph the tested spiral bevel pinion. 

Many of his photographs are contained in this volume. 


Available from 


NASA Center for Aerospace Information 
7121 Standard Drive 
Hanover, MD 21076 
Price Code: A06 


National Technical Information Ser\dce 
5285 Port Royal Road 
Springfield, VA 22100 
Price Code: A06 



TABLE OF CONTENTS 


CHAPTER ONE: INTRODUCTION 1 

1.1 Background 1 

1 .2 Numerical Analyses of Gears 3 

1 .3 Overview of Chapters 5 

CHAPTER TWO: GEAR GEOMETRY AND LOADING 7 

2.1 Introduction 7 

2.2 Basics of Spiral Bevel Gear Geometry 7 

2.3 Teeth Contact and Loading of a Gear Tooth 1 1 

2.4 Gear Materials 16 

2.5 Motivation to Model Gear Failures 16 

2.5. 1 Gear Failures 1 8 

2.5.2 OH-58 Spiral Bevel Gear Design Objectives 19 

2.6 Chapter Summary 1 9 

CHAPTER THREE: COMPUTATIONAL FRACTURE MECHANICS 21 

3.1 Introduction 21 

3.2 Fracture Mechanics and Fatigue 2 1 

3.2.1 Fatigue 23 

3.2.2 Example: Two dimensional, mode I dominant fatigue crack growth 

simulation with static, proportional loading 27 

3.2.3 Example: Three dimensional, mode I dominant fatigue crack growth 

simulation with static, proportional loading 3 1 

3.3 Fracture Mechanics Software 33 

3.4 Chapter Summary 34 

CHAPTER FOUR: FATIGUE CRACK GROWTH RATES 35 

4.1 Introduction 35 

4.2 Fatigue Crack Closure Concept 35 

4.3 Application of Newman’ s Model to AISI 93 1 0 Steel 40 

4.4 Sensitivity of Growth Rate to Low R 44 

4.6 Chapter Summary 46 

CHAPTER FIVE: PREDICTING FATIGUE CRACK GROWTH 
TRAJECTORIES IN THREE DIMENSIONS UNDER MOVING, NON- 
PROPORTIONAL LOADS 47 

5.1 Introduction 47 

5.2 BEM Model 47 

5.2.1 Loading Simplifications 49 

5.2.2 Influence of Model Size on SIF Accuracy 5 1 

5.3 Initial SEF History Under Moving Load 54 


NAS A/CR— 2000-2 1 0062 


V 



5.4 Method for Three Dimensional Fatigue Crack Growth Predictions Under Non- 

Proportional Loading 58 

5.4. 1 Literature Review 58 

5.4.2 Proposed Method 59 

5.4.3 Approximations of Method 63 

5.5 Simulation Results 64 

5.6 Chapter Summary 67 

CHAPTER SIX: EXPERIMENTAL RESULTS 69 

6.1 Introduction 69 

6.2 Test Results 69 

6.3 Fractography 71 

6.3.1 Overview 71 

6.3.2 Results 73 

6.4 Chapter Summary 79 

CHAPTER SEVEN: DISCUSSION AND SENSITIVITY STUDIES 81 

7.1 Introduction 81 

7.2 Comparisons of Crack Growth Results 81 

7.3 Sensitivity Studies 85 

7.3.1 Fatigue Crack Growth Rate Model Parameters 86 

7.3.2 Crack Closure Model Parameters 87 

7.3.3 Loading Assumptions 89 

7.4 Highest Point of Single Tooth Contact (HPSTC) Analysis 96 

7.5 Chapter Summary 99 

CHAPTER EIGHT: CONCLUDING REMARKS 101 

8.1 Accomplishments and Significance of Thesis 101 

8.2 Recommendations for Future Research 103 

APPENDIX A 104 

APPENDIX B 106 

APPENDIX C 108 

REFERENCES 110 


NAS A/CR— 2000-2 10062 


VI 



LIST OF ABBREVIATIONS 


AGMA 

American Gear Manufacturers Association 

BEM 

Boundary element method 

EDM 

Electro-discharge machined 

FEM 

Finite element method 

FRANC3D 

FRacture ANalysis Code - 3D 

HPSTC 

Highest point of single tooth contact 

LEFM 

Linear elastic fracture mechanics 

NASA/GRC 

National Aeronautics and Space Administration - Glenn Research 
Center 

OSM 

Object Solid Modeler 

RC 

Rockwell C 

SEM 

Scanning electron microscope 

SIF 

Stress intensity factor 


vii 


NAS A/CR— 2000-2 1 0062 




CHAPTER ONE: 
INTRODUCTION 


1.1 Background 

A desirable objective in the design of aircraft components is to minimize the 
weight. A lighter aircraft operates more efficiently. A helicopter’s transmission 
system is one example where design is focused on weight minimization. A 
transmission system utilizes various types of gears, such as spur gears and spiral bevel 
gears. Because spur gear geometry is relatively simple, optimizing the design of these 
gears using numerical methods has been researched significantly. However, the 
geometry of spiral bevel gears is much more complex, and less research has focused 
on using numerical methods to evaluate their design and safety. 

One obvious method to minimize the weight of a gear is to reduce the amount 
of material. However, removing material can sacrifice the strength of the gear. In 
addition, fatigue cracks in gears are a design concern because of the cyclical loading 
on a gear tooth. Research shows that the size of a spur gear’s rim with respect to its 
tooth height determines the crack trajectories [Lewicki et al. 1997a, 1997b]. This 
knowledge is critical because it allows the designer to predict failure modes based on 
geometry. 

Two common failure modes of a gear are rim fracture and tooth fracture. Rim 
fracture, shown in Figure 1.1 [Albrecht 1988], can be catastrophic and lead to the loss 
of the aircraft and lives. On the other hand. Figure 1.2 is an example of a tooth 
fracture [Alban 1985]. Tooth fracture is the benign failure mode because it is most 
often detected prior to catastrophic failure. Knowing how crack trajectories are 
affected by design changes is important with respect to these two failure modes. 



Figure 1.1: Spiral bevel gear rim failure [Albrecht 1988]. 
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Figure 1.2: Spiral bevel gear tooth failure [Alban 1985], 

In general, gears in rotorcraft applications are designed for infinite life; the 
gears are designed to prevent any type of failure from occurring. Developing a 
damage tolerant design approach could reduce cost and increase effectiveness of the 
gear. Lewicki et al.’s work on determining the effect of gear rim thickness on crack 
trajectories is a good example of how damage tolerance can be applied to gears. 
Knowing how the gear’s geometry affects the failure mode allows a designer to select 
a geometry such that, if a crack were to develop, the failure mode would be benign. 
Other examples of damage tolerant design can be found in aircraft structures [Swift 
1984] [Rudd 1984] [Miller et al. 1999], helicopter rotor heads [Irving et al. 1999], and 
train rails [Jeong et al. 1997]. 

Damage tolerance involves designing under the assumption that flaws exist in 
the structure [Rudd 1984]. The initial design then focuses on making the structure 
sufficiently tolerant to the flaws such that the structural integrity is not lost. Damage 
tolerant design allows for multiple load paths to prevent the structure from failing 
within a specified time after one element fails. In this regard, gears would be designed 
for the benign failure mode, tooth failure, as opposed to rim failure, which could be 
catastrophic. 

Current American Gear Manufacturers Association (AGMA) standards use 
tables and indices to approximate the strength characteristics of gears [AGMA 1996]. 
The finite element method (FEM) and boundary element method (BEM) are becoming 
more useful and common approaches to study gear designs. A primary reason for this 
is the tremendous increase in computing power. Section 1.2 summarizes recent 
research related to modeling gears numerically. 

Limited work has focused on predicting crack trajectories in spiral bevel gears. 
This is most likely because a spiral bevel gear’s geometry is complex and requires a 
three dimensional representation. Structures with uncomplicated geometries, such as 
spur gears, can be modeled in two dimensions. Modeling an object in three 

dimensions requires a crack to also be modeled in three dimensions. Three 
dimensional crack representations introduce unique challenges that do not arise when 
modeling in two dimensions. 

A three dimensional crack model consists of a continuous crack front. When a 
simpler geometry allows for a two dimensional simplification, a crack front is now 
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represented by a single point, the crack tip. At a crack tip there are only two modes of 
displacement; in three dimensional models, however, there is a distribution of three 
modes of displacement along the crack front. Propagating a crack in two dimensions 
is completely defined by a single angle and extension length. On the other hand, along 
the crack front there is a distribution of angles and lengths. 

Codes developed by the Cornell Fracture Group at Cornell University, such as 
Object Solid Modeler (OSM) and FRacture ANalysis Code - 3D (FRANC3D), have 
been developed to handle three dimensional fracture problems. FRANC3D explicitly 
models cracks and predicts crack trajectories under static loads. The crack growth 
models are based on accepted fatigue crack growth and linear elastic fracture 
mechanics (LEFM) mixed mode theories. 

Because gears operate at high loading frequencies, the actual time from crack 
initiation to failure is limited. As a result, crack trajectories and preventing 
catastrophic failure modes are the primary concern in gear design. Crack growth rates 
are not as important. The goal of this research is to investigate issues related to 
predicting three dimensional fatigue crack growth in spiral bevel gears. A simulation 
that allows for arbitrarily shaped curved crack fronts and crack trajectories will be 
most accurate. In addition, the loading on a tooth as a function of time, position, and 
magnitude should be considered. 

1.2 Numerical Analyses of Gears 

Computational fracture mechanics applied to gear design is a relatively novel 
research area. As a result, the majority of work has been limited to two dimensional 
analyses. In three dimensions, very little work has predicted crack trajectories in 
gears. This section summarizes some pertinent developments in applying numerical 
methods and fracture mechanics to gear design. 

The complexity of two dimensional gear analyses has evolved. Albrecht 
[1988] used the FEM to investigate gear tooth stresses, gear resonance, and 
transmission noise. Individual gear teeth were modeled in two dimensions and the 
increase in accuracy when using the FEM over AGMA standard indices for 
calculating gear tooth root stresses was demonstrated. Blarasin et al. [1997] used the 
FEM and weight function technique to evaluate stress intensity factors (SIFs) in 
specimens similar to spur gear teeth. Cracks with varying depths were introduced in 
two dimensional models and a constant single point load was applied. The SIFs were 
determined as a function of crack depth. Fatigue lives were calculated, but predictions 
of the crack trajectory were never performed. Flasker et al. [1993] used two 
dimensional FEM to analyze fatigue crack growth in a gear of a car gearbox. The 
analyses considered highest point of single tooth contact (HPSTC), but variable 
loading at that point. Residual stresses from the case and core were simulated with 
thermal loading. Based on a given load history, the crack was incrementally 
propagated. Lewicki et al. [1997a, 1997b] combined FEM and LEFM to investigate 
crack trajectories in thin rimmed spur gears. The work successfully matched crack 
trajectory predictions to experiments. 

Limited three dimensional crack analyses of gears have been achieved. The 
work most often concerns simple geometries and loading conditions. Pehan et al. 
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[1997] used the FEM to look at two and three dimensional spur gear models. Residual 
stresses due to case hardening were modeled as nodal thermal loads. Two different 
sized models were analyzed: one tooth including the arc length of the gear rim directly 
below the tooth and three teeth with the corresponding gear rim arc length. To 
determine the new crack front, they used a criterion such that the SIFs along the new 
front should be constant. Paris’ model was used to calculate the fatigue lives based on 
the SIFs near the midpoint of the crack front. A constant load location with constant 
magnitude and simple spur gear geometry allowed Pehan et al. to consider only crack 
opening (mode I) effects. Their method for determining the new crack fronts is 
computationally intensive and limited since three dimensional effects are not 
accommodated. 

Lewicki et al. [1998] performed three dimensional crack propagation studies 
using the FEM and BEM to investigate fracture characteristics of a split tooth gear 
configuration. The geometry of the split tooth configuration is similar to a spur gear. 
The analyses used single load locations and explored propagation paths for various 
crack locations. The strong point of this work is that three dimensional simulations of 
crack trajectories were performed in addition to calculating fatigue crack growth rates. 

Very little work, in addition to Lewicki et al.'s research, has used the BEM to 
analyze gears. Sfakiotakis [1997] performed two dimensional BEM analyses of gear 
teeth considering mechanical and thermal loads. Rather than perform trajectory 
predictions, they calculated SIFs for different size initial cracks with various loading 
conditions and crack locations. Fatigue loading was not considered. Fu et al. [1995] 
also used the BEM for stress analysis related to optimizing the forging die of spiral 
bevel gears. 

The progression of research related to computer analysis of gears has led to the 
investigation of crack growth in spiral bevel gears. FEM models of spiral bevel gears 
can be created from Handschuh et a/. ’s [1991] computer program that models the 
cutting process of spiral bevel gears to determine tooth surface coordinates in three 
dimensions. Litvin et al. [1996] utilized this program, in conjunction with tooth 
contact analysis [Litvin et al. 1991], to determine how bearing (contact between 
mating gear teeth) changes with different spiral bevel gear tooth surface designs. 
Transmission error curves were generated that gave an indication of the efficiency of 
the gear, 

Along with Litvin et al.'s work [1991], tooth contact analysis of mating gears 
has been explored by Bibel et al. [1995 and 1996], Savage et al. [1989], and Bingyuan 
et al. [1991]. Bibel et al. successfully modeled multi-tooth spiral bevel gears with 
deformable contact using the FEM. They conducted a stress analysis of mating spiral 
bevel gears and analytically modeled, using gap elements from general purpose finite 
element codes, the rolling contact between the gear teeth. Bibel et al.'s work can be 
used to investigate how changes in gear geometry affect tooth deflections. Variations 
in tooth deflections can alter the contact zone between gear teeth. Savage et al. 
developed analytical methods to predict, using tooth contact analysis, the shift in 
contact ellipses due to elastic deflections of a spiral bevel gear’s shafts and bearings 
under loads. Savage et al. and Bibel et al.'s work was related to spiral bevel gears, 
however, they did not incorporate fracture mechanics. On the other hand, Bingyuan et 
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al. approximated the geometry of gears in contact as a pair of disk rollers compressed 
together. The linear elastic stresses in the disks could be written in closed form. The 
SIFs were calculated using the closed form expressions. Bingyuan et aL's primary 
focus was to calculate surface fatigue life and compute crack growth rates. No 
trajectory predictions were made. 

TTie majority of the aforementioned research on spiral bevel gears is unrelated 
to failure, but rather associated with design and efficiency; methods have been 
developed to create numerical models of spiral bevel gears and predict contact areas. 
Crack trajectories have been predicted in gears with simpler geometry that can be 
represented by two dimensional models. This thesis is a natural extension of the 
research to date. The next step is to computationally model fatigue crack trajectories 
in spiral bevel gears. 

1.3 Overview of Chapters 

This thesis is divided into eight chapters. The first and last chapters are 
overview and summary. The remaining chapters each build upon one another and 
propiose, apply, and evaluate methods for predicting fatigue crack growth in spiral 
bevel gears. 

Chapter Two contains background information on gears, with particular 
attention to spiral bevel gears. The objective is to define vocabulary and concepts 
related to spiral bevel gears that will be used throughout the thesis. In addition, the 
work of the thesis is further motivated by examples of gear failures and the current 
design objectives for gears. 

A focus of this thesis is to demonstrate that computational fracture mechanics 
can be used to analyze complex gear geometries under realistic loading conditions. 
LEFM and fatigue theories that are utilized to accomplish this task are presented in 
Chapter Three. Methods that are currently implemented in two and three dimensions 
to compute crack trajectories are demonstrated through examples. 

Chapter Four explores the significance of compression loading on calculated 
crack growth rates. The magnitude of compressive stresses in a gear’s tooth root is a 
function of the rim thickness. If fatigue crack growth rates are highly sensitive to this 
compression, then growth rates may warrant more attention in designing gears. The 
concept of fatigue crack closure is used to investigate fatigue crack propagation rates 
in AISI 9310, a common gear steel. First, the concept of fatigue crack closure is 
discussed. A material-independent method is presented for obtaining fatigue crack 
growth rate data that do not vary with stress ratio. The method is demonstrated using 
data at various stress ratios for pressure vessel steel. Next, the concepts are applied to 
AISI 9310 steel data to obtain an intrinsic fatigue crack growth model. This model is 
used to investigate the effect of low stress ratios on fatigue crack growth in AISI 9310. 

Chapter Five is an initial investigation into predicting three dimensional 
fatigue crack trajectories in a spiral bevel pinion under a moving load. First, a 
boundary element model of a pinion is developed. A method to represent the moving 
contact area on a gear tooth is discussed. Next, studies are conducted to determine the 
smallest model that still accurately represents the operating conditions of the pinion. 
Once the model is defined, a crack is introduced into the model, and the initial stress 
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intensity factor history under the moving load is calculated. A method to predict 
fatigue crack trajectories under the moving load is proposed. The method is then 
applied to predict fatigue crack growth trajectories and rates in a spiral bevel pinion. 

Fatigue crack growth results from a spiral bevel pinion in operation are 
necessary to validate the predictions. The sponsor of the research efforts of this thesis, 
NASA-Glenn Research Center (NASA/GRC), provided a pinion that was tested in 
their gear test fixture. Notches were fabricated into several of the teeth’s roots prior to 
beginning the test. The test data and crack growth results are presented in Chapter 6. 
In addition, in an effort to obtain crack front shape and crack growth rate information, 
the fracture surfaces are obser\'ed with a scanning electron microscope, and the results 
are given in the chapter. 

The crack trajectory and fatigue life results from the simulation and the tested 
pinion are compared in Chapter Seven. To gain insight into the discrepancies between 
the prediction and test, the influence of model parameter assumptions and loading 
simplifications on crack trajectories and calculated fatigue crack growth rates are 
studied. Next, the necessity of the moving, non-proportional load crack growth 
method is evaluated by comparing the results to predictions that assume proportional 
loading. 

Finally, Chapter Eight summarizes the accomplishments of the work in the 
previous chapters. Implications of the research conducted and suggestions for future 
work are given. 
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CHAPTER TWO: 

GEAR GEOMETRY AND MODELING 


2.1 Introduction 

Chapter Two covers the basic terms and geometry aspects of a spiral bevel gear. 
This terminology and background is essential to motivate the numerical simulations of 
this thesis. A gear’s design and geometry can be quite complex; however, only the 
fundamentals are explained in this chapter 

2.2 Basics of Spiral Bevel Gear Geometry 

Gears are used in machinery to transmit motion. Gears operate in pairs. The 
two mating gears have similar shapes. The smaller of the mating gears is called the 
pinion, and the larger the gear. Motion is transferred from one gear to another by 
successively engaging teeth. 

There are various types of gears. The shape of the teeth and the angle at which 
the mating gears are mounted are a few of the distinguishing characteristics between 
the gear types. Gears with intersecting shafts are called bevel gears. The most 
common angle to mount bevel gears is 0 = 90° , although any intersecting angle could 
be used. A bevel gear’s form is conical. For comparison, as illustrated in Figure 2.1, 
spur gears are cylindrical, and the shafts of the gears are parallel. The geometry of a 
spur gear can be almost fully illustrated in two dimensions. However, the conical 
shape of a bevel gear requires a three dimensional illustration. This two and three 
dimensional difference is where the complexity of the work contained in this thesis 
lies. 



a) Spur gears operate with parallel axes 
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B 

b) Bevel gears operate with intersecting axes 

Figure 2. 1 ; Schematics of spur (a) and bevel (b) gears. 

The cone defined by the angle between a bevel gear’s axis and the line of 
tangency with the mating gear is called the pitch cone. In Figure 2.1b, 6i and ft 
define the pitch cones. The gear ratio is the ratio of the angular frequencies of the 
mating gears, 0hl(O\, which also equals the ratio of sin(ft) to sin(ft), or, due to 
geometry, the ratio of the number of gear teeth to the number of pinion teeth. 




Figure 2.; 


b) Spiral bevel gear 

me 


Two common bevel gears are the straight bevel gear and the spiral bevel gear. 
The main difference between these two gears is the shape of their teeth. The teeth of 
the straight bevel gear are straight, and the teeth of the spiral bevel gear are curved. 
Figure 2.2 illustrates this difference. When looking along the axis of a spiral bevel 
gear, the teeth will either curve counterclockwise or clockwise, depending on whether 
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the gear is left- or right-handed, respectively. So that the teeth can fit together, or 
mesh, a spiral bevel gear and pinion will always have opposite hands. The thickness 
and height of a spiral bevel gear tooth varies along the cone. The larger end of the 
tooth is the heel, and the smaller the toe. The curvature of the tooth creates concave 
and convex tooth surfaces on opposite sides of the tooth. Figure 2.3. 



The tooth profile, as shown in Figure 2.4, is one side of the cross section of a 
gear tooth. The fillet curve is at the bottom of the tooth profile where it joins the space 
between the teeth. The region of the tooth near the fillet is the bottom land, and the 
area near the top of the profile is the top land. 


Top Land 



Figure 2.4: Schematic of cross section of a gear tooth. 

The advantage of the spiral bevel gear’s curved teeth is to allow for more than 
one tooth to be in contact at a time. This makes it significantly stronger than a straight 
bevel gear of equal size. Consequently, spiral bevel gears are commonly found in 
high speed and high force applications. One such application, which is the focus of 
this thesis, is in helicopter transmission systems. The mating spiral bevel gears in the 
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transmission system convert the power from the horizontal engine shaft to the vertical 
shaft of the main rotor. Gears in this application typically operate at rotational speeds 
of 6000 rpm and transmit on the order of 300 hp of power. 

Many parallel axis gears, such as spur gears, have involute tooth profiles. As 
sketched in Figure 2.5, the involute curve can be visualized by unwrapping thread 
from a spool while keeping the thread taut. The path traced by the end of the string is 
an involute curve. The spool is the evolute curve. All involute gear geometries are 
generated from circle evolute curves. The involute curve then becomes a spur gear 
tooth’s profile. A closed form solution for the coordinates along the curve exists for 
this type of geometry. As a result, the tooth’s surface coordinates can be calculated 
with relative ease. However, the geometry of a spiral bevel tooth is much more 
complex, and there is no closed form solution to describe the surface coordinates. 
Handschuh et al. [1991] developed a program to numerically calculate the surface 
coordinates of a spiral bevel gear tooth. The program models the kinetics of the 
cutting process in creating the gear, along with the basic gear geometry. The program 
calculates the coordinates of a spiral bevel gear tooth in three dimensions for use as 
input to a finite element model. The numerical models in this thesis were all created 
using the tooth geometry coordinates as defined by Handschuh et o/.’ s program. 


Involute Curv^e 



Figure 2.5: Generation of an involute curve. 
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Figure 2.6: OH-58 spiral bevel pinion with two fractured teeth. 


A spiral bevel gear set is used in the main rotor transmission of the U.S. 
Army’s OH-58 Kiowa Helicopter. An OH-58 spiral bevel pinion that exhibited tooth 
fracture during an experiment is shown in Figure 2.6. 

The geometry of the OH-58 gear set will be used throughout this thesis. In the 
set, a 19 tooth spiral bevel pinion meshes with a 71 tooth spiral bevel gear. The 
pinion’s shafts are supported by ball bearings. The input torque is applied at the end 
of the pinion’s large shaft. The approximate dimensions of a pinion tooth are given 
schematically in Figure 2.7. 



Figure 2.7: Approximate dimensions of OH-58 spiral bevel pinion tooth. 


2.3 Teeth Contact and Loading of a Gear Tooth 

According to the theory of gears, there is a point of contact between a spiral 
bevel gear and pinion at any instant in time where their surfaces share a common 
normal vector. In reality, the tooth surfaces deform elastically under the contact. The 
deformation spreads the point of contact over a larger area. The larger area has 
traditionally been approximated using Hertzian contact theory. This contact is 
conventionally idealized to spread over an elliptical area [Johnson 1985]. The center 
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of the ellipse is the mean contact point, which determines the contact ellipse’s location 
on the tooth surface. The orientations of the ellipse’s minor and major axes are 
defined by the tooth surface’s geometry, curvature, and the alignment between the 
gear and pinion. The length of the axes is a function of the load. It can be shown that 
the ratio of the axes’ lengths is constant and is not _a function of the load. The form of 
the equations for the length of the ellipse’s semi-major and semi-minor axes, a and h, 
respectively, is [Johnson 1985] [Timoshenko cr a/. 1970]: 


« = / 


3;r 


1 / 


(2.1a) 


b^g 


37T 


-. 1 / 


(2.1b) 


where / and g are functions defined by the geometry. The magnitude of force, F, 
exerted on the tooth is proportional to the input torque level and gear geometry. 

The meshing of the mating gear teeth is a continuous process. The position of 
the area of contact and magnitude of the force exerted between the teeth varies with 
time as the gear rotates. Figure 2.8 illustrates schematically the progression of the 
contact area along a tooth of a left-handed spiral bevel pinion. In the schematic, the 
continuous process has been discretized into a series of elliptical contact patches, or 
load step increments. The darkened arrow demonstrates the direction the load moves. 
The actual tooth contact pattern during operation is a function of the alignment of the 
gear and pinion. 



left-handed spiral bevel pinion tooth. 


Overlap in tooth contact between adjacent teeth results in two modes of contact: 
single tooth contact and double tooth contact. At the beginning of a meshing cycle for 
one tooth, two teeth of the pinion are in contact with the gear. As the pinion rotates, 
the adjoining tooth loses contact_with the gear and only one pinion tooth receives all of 
the force. As the pinion continues to rotate, the load moves further up the pinion 
tooth, and the next pinion tooth comes into contact with the gear; the force on a pinion 
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tooth is once again reduced due to the double tooth contact. The contact area will 
differ for single tooth and double tooth contact. The change in area of the contact is 
schematically illustrated in Figures 2.8 and 2.9. 


Time Step Tooth 1 Tooth 2 



Figure 2.9: Schematic of load progression on adjacent pinion teeth. 
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In Figure 2.9, tooth 1 and 2 are two adjacent teeth of a spiral bevel pinion. The 
ellipses represent “snap shot” areas of contact between a gear and a pinion’s tooth. 
The darkened ellipse is the area that is currently in contact with the gear at a particular 
instant in time. Similar to Figure 2.8, the larger ellipses represent single tooth contact, 
and the smaller are areas of double tooth contact. The first row in Figure 2.9 begins 
with tooth 1 at the last moment of single tooth contact. After a discrete time step, the 
load on tooth 1 has progressed up the tooth and tooth 2 has come into contact near the 
root, as depicted in row two. In the final row, or time step, tooth 1 loses contact and 
tooth 2 advances into the stage of single tooth contact. 

It is seen in Figures 2.8 and 2.9 that the contact area between mating spiral 
bevel gear teeth moves in three spatial dimensions during one load cycle. Most of the 
previous research into numerically calculating crack trajectories in gears has been 
performed on spur gears with two dimensional analyses and has not incorporated the 
moving load discussed above. Instead, a single load location on the spur gear tooth 
that produces the maximum stresses in the tooth root during the load cycle has been 
used to analyze the gear. This load position corresponds to the highest point of single 
tooth contact (HPSTC). Contact between spur gear teeth only moves in two 
directions, and, therefore, this simplification to investigate a spur gear under a fixed 
load at the HPSTC has proven successful [Lewicki 1995] [Lewicki et al. 1997a]. 
However, since the contact area between mating spiral bevel gear teeth moves in three 
dimensions, the crack front trajectories could be significantly influenced by this three 
dimensional effect. As a result, trajectories under the moving load should be predicted 
first and compared to trajectories considering only a fixed loading location at HPSTC. 
This approach is detailed in Chapters 5 and 7. 

It has been implicitly assumed in the above discussion that the traction, or 
force over the contact area, is normal to the surface. Dike [1978] points out that this 
assumption is valid if there are no frictional forces in the contact area. He also states 
that is the case with gears since a lubricant is always used. The lubricant will make 
the magnitude of the frictional forces small competed to the normal forces. This 
assumption will be utilized in the numerical simulations. 

In the same paper. Dike also asserts that there are two main areas in a gear 
tooth where the bending stresses may cause damage. The first is the location of 
maximum tensile stresses at the fillet of the tooth on the same side as the load. The 
second is at the fillet of the tooth on the side opposite the load, where the maximum 
compressive stresses occur. 

This can be visualized by drawing an analogy between a cantilever beam and a 
gear tooth. Figure 2.10, Basic beam theory predicts that the maximum tensile stress 
occurs at the beam/wall connection on the outer most fibers on the same side as the 
applied load. The maximum compressive stress occurs at the same vertical location, 
on the side opposite the load. Similarly, as a gear tooth is loaded, it creates tensile 
stresses in the tooth root of the loaded side. In the root of the side opposite the load, 
there are compressive stresses. These compressive stresses might also extend into the 
fillet and root of the next tooth. 
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Maximum 

compressive 

stress 



Maximum 

compressive 

stress 



Figure 2.10: Stresses in cantilever beam (a) are analogous to gear tooth root (b). 

The compressive stresses are noteworthy because Lewicki et al. [1997b] 
showed that the magnitude of the compressive stress increases as a gear’s rim thickness 
decreases. The compressive stress could affect the crack propagation trajectories and 
crack growth rates. However, it is demonstrated in Chapter 4 that low stress ratios, i.e. 
large compressive stresses compared to tensile stresses, do not have a significant 
influence on crack predictions. 

Up to this point, only frictional loads and traction normal to the tooth’s surface 
have been discussed. The normal loads are the only loading conditions to be 
considered in this thesis. However, additional sources do produce forces on the gear. 
Some of these additional loads include dynamic effects, centrifugal forces, and 
residual stresses due to the case hardening of the gear. In addition, since a lubricant is 
always used when gears are in operation, lubricant could get inside a crack and create 
hydraulic pressure. 
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2.4 Gear Materials 

As discussed in Section 2.2, spiral bevel gears are commonly used in helicopter 
transmission systems. In this application, the gear’s material impacts the life and 
performance of the gear. Most often a high hardenable iron or steel alloy is used. The 
traditional material for the OH-58 spiral bevel gear is AISI 9310 steel (AMS 6265 or 
AMS 6260). Some other aircraft quality gear steels are VASCO X-2, CBS 600, CBS 
1000, Super Nitroalloy, and EX-53. The choice of material is dependent on operating 
variables such as temperature, loads, lubricant, and cost. The material characteristics 
most important for gears are surface fatigue life, hardenability, fracture toughness, and 
yield strength. Table 2.1 shows the chemical composition of AISI 9310 [AMS 1996]. 
Table 2.2 contains relevant material properties. 


Table 2 . 1 : Chemical composition of AISI 93 10 by weigl 



C 

Mn 

P 

S 

Si 

Cu 

Ni 

Cr 

B 

Mo 

Fe 

Minimum 

0.07 

0.40 

-- 

— 

0.15 

— 

3.00 

1.00 

— 

0.08 

95.30 

Maximum 

0.13 

0.70 

0.015 

0.015 

0.35 

0.35 

3.50 

1.40 

0.001 

0.15 

93.39 


ht percent [AMS 1996]. 


Most gears are case hardened. Case hardening increases the wear life of the gear. 
In general, the gears are vacuum carburized to an effective case depth’ of 0.032 in - 
0.040 in (0.813 mm - 1.016 mm). The case hardness specification is 60 - 63 Rockwell 
C (RC), and the core hardness is 31 - 41 RC [AGMA 1983]. 


Table 2.2: Material properties of AISI 9310. 


^ 

Tensile Strength" 

185 X psi 

Yield Strength' 

160 X 10"' psi 

Young’s Modulus 

30 X 10^ psi 

Poisson’s Ratio 

0.3 

Fracture Toughness^ 

85 ksi*in" 

Average Grain Size'* 

ASTM No. 5 or finer 
(< 0.00244 in) 


2.5 Motivation to Model Gear Failures 

Gear failures can be categorized into several failure modes. Tooth bending, 
pitting, spalling, and thermal fatigue can all be placed in the category of fatigue 
failures. Examples of impact type of failures are tooth shear, tooth chipping, and case 
crushing. Wear and stress rupture are two additional modes of failure. According to 
[Dudley 1986], the three most common failures are tooth bending fatigue, tooth 
bending impact, and abrasive tooth wear. He gives examples of a variety of failures 
from tooth bending fatigue to spalling to rolling contact fatigue in both spur and spiral 
bevel gears. 


' The effective case depth is defined as the depth to reach 50 RC. 
‘ [Coy eta/. 1995] 

■ [Townsend e/ a/. 1991] 

^ [AMS 19961 
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The focus of this thesis is on tooth bending fatigue failure because this is one 
of the most common failures. In general, tooth bending fatigue crack growth can lead 
to two types of failures. In rotorcraft applications, the type of failure could be either 
benign or catastrophic. Crack propagation that leads to the loss of one or more 
individual teeth will most likely be a benign type of failure. The remaining gear teeth 
will still be able to sustain load, and the failure should be detected due to excessive 
vibration and noise. On the other hand, a crack that propagates into and through the 
rim of the gear leaves the gear inoperable. The gear will no longer be able to carry 
any load, and will most likely lead to loss of aircraft and life. 

Alban [1985, 1986] proposes a “classic tooth-bending fatigue” scenario. He 
suggests five conditions that characterize the “classic” failure: 

1 . The origin of the fracture is on the concave side in the root. 

2. The origin is midway between the heel and the toe. 

3. The crack propagates first slowly toward the zero-stress point in the root. 
As the crack grows, the location of the zero-stress point moves toward a 
point under the root of the convex side. The crack then progresses outward 
through the remaining ligament toward the convex side’s root. 

4. As the crack propagates, the tooth deflection increases only up to a point 
when the deflection is large enough that the load is picked up 
simultaneously by the next tooth. Since the load on the first tooth is 
relieved, the rate of increase in the crack growth rate decreases. 

5. No material flaws are present. 

Alban presents results from a photoelastic study of mating spur gear teeth. The 
study demonstrates the shift in the zero-stress point. The zero-stress point is where the 
tensile stresses in the root of loaded side of the tooth shift to compressive stresses on 
the load free side. Figure 2.11 shows stress contours for two mating spur gear teeth. 
In the bottom gear, one of the teeth is cracked and another tooth has already fractured 
off. The teeth of the top gear are not flawed. By comparing contours between the 
mating cracked and uncracked teeth, it is easy to pick out the zero-stress location shift 
toward the root of the load free side. The shift of the zero-stress location demonstrates 
the changing stress state in the tooth. This changing stress state drives the crack to 
turn. The point in the two dimensional cross section where the crack turns is actually 
a ridge when the third spatial dimension, the length of the tooth, is considered. This 
classic tooth failure scenario will be used as a guideline when evaluating the 
prediction and experimental results in the following chapters. 
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Zero-stress point „ ^ 

Compressive stress ^ Tensile stress 



Crack 


Figure 2.11: Photoelastic results from mating spur gear teeth (stress contour 
photograph from [Alban 1985]). 


2.5.1 Gear Failures 

Gears in rotorcraft applications are currently designed for infinite life. 
Therefore, gear failures are not common. However, failures do occur primarily as a 
result from manufacturing flaws, metallurgical flaws, and misalignment. 

Dudley [1996] gives an overview of the various factors affecting a gear’s life. 
Some of the more common metallurgical flaws listed are case depth too thin or too 
thick, grinding bums on the case, core hardness too low, inhomogeneities in the 
material microstmcture, composition of the steel not within specification limits, and 
quenching cracks. In addition, examples of surface durability problems, such as 
pitting, are presented. A pitting flaw could develop into a starter crack for a fatigue 
failure. riL: 

Pepi [1996] examined a failed spiral bevel gear in an Army cargo helicopter. 
A grinding bum was determined as the origin of the fatigue crack. In addition, it was 
learned that the carburized case was deeper than acceptable limits in the area of the 
crack origin, which contributed to crack growth. Roth et al. [1992] determined a 
microstructure inhomogeneity, introduced during the remelting process, to be the 
cause of a fatigue crack in a carburized AISI 9310 spiral bevel gear. Both of these 
failures could be classified as manufacturing flaws. 

Albrecht [1988] gives an example of a series of failures in the Boeing Chinook 
helicopter, which were caused by gear resonance with insufficient damping. Couchon 
et al. [1993] gives an example of a gear failure resulting from excessive misalignment. 
The excessive misalignment was due to a failed bearing that supported the pinion. 
The misalignment led to a fatigue crack on the loaded side of the tooth. An analysis of 
an input spiral bevel pinion fatigue crack failure in a Royal Australian Navy helicopter 
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is given by McFadden [1985]. These examples demonstrate that gear failures do 
occur in service. 

Gear experts are researching ways to make gears quieter and lighter through 
changes in the geometry. However, at the same time there is a tradeoff between 
weight, noise, and reliability. Geometry changes could have negative effects on the 
strength and crack trajectory characteristics of the gear. A design tool to predict the 
performance of proposed gear designs and changes, such as discussed by Lewicki 
[1995], would be extremely useful. Savage er al. [1992] used an optimization 
procedure to design spiral bevel gears using gear tooth bending strength and contact 
parameters as constraints. Including effects of geometry changes on the strength and 
failure modes could contribute greatly to his procedures. 

2.5.2 OH-58 Spiral Bevel Gear Design Objectives 

In rotorcraft applications, a primary source of vibration of the gear box is 
produced by the spiral bevel gears [Coy et al. 1987] [Lewicki ei al. 1993]. In turn, the 
vibration of the gear box accounts for the majority of the interior cabin noise. As a 
result, recent design has focused on modifying the gear’s geometry to reduce the 
vibration and noise. In addition, due to the application of the gear, a continuous 
design objective is to make the gear lighter and more reliable. 

Adjusting the geometry of the gear, however, may jeopardize the gear’s 
strength characteristics. Lewicki et al. [1997a] showed that the failure mode in spur 
gears is closely related to the gear’s rim thickness. It was demonstrated that if an 
initial flaw exists in the root of a tooth, the crack would propagate either through the 
rim or through the tooth for a thin rimmed and thick rimmed gear respectively. As a 
result, a tool to evaluate the strength and fatigue life characteristics of proposed gear 
designs would be useful. 

Albrecht [1988] demonstrated that AGMA standards to determine gear stresses 
and life were insufficient. He also showed the advantages of a numerical simulation 
method, such as the FEM, over the currently accepted AGMA standards at that time. 
The work of this thesis is an extension of the numerical approaches to determine gear 
stresses and life. 

2.6 Chapter Summary 

This chapter covered basic terminology and geometry aspects of gears. 
Concepts related to spiral bevel gears were the primary focus. In addition, methods to 
visualize and model the contact between mating spiral bevel gears were presented. 
Characteristics of a common gear steel, AISI 9310, were summarized. These 
materials properties will be used in the numerical simulations. Finally, some 
examples of gear failures and gear design objectives were discussed to motivate the 
significance of modeling gear failures numerically. 
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CHAPTER THREE: 

COMPUTATIONAL FRACTURE MECHANICS 


3.1 Introduction 

This chapter discusses areas of computational fracture mechanics relevant to 
the work of this thesis. The areas of focus are LEFM, fatigue, and the BEM. The 
BEM is used in a fashion similar to the more common FEM. The primary difference 
between the methods in three dimensional elasticity problems is that with the BEM 
only the surfaces, or boundaries, are meshed, as opposed to the volume that is meshed 
in the FEM. In computational LEFM, the displacement and/or stress results from a 
numerical analysis are used to calculate the SIFs. The SIFs are in turn used to predict 
how and where a crack may grow. 

The analyses of this work are conducted using a suite of computational fracture 
mechanics programs developed by the Cornell Fracture Group. OSM is used to create 
a geometry model of the OH-58 spiral bevel pinion. FRANC3D is used as a pre- and 
post-processor to the boundary element solver program, BES. FRANC3D has built in 
features to compute SIFs using the displacement correlation technique. 

3.2 Fracture Mechanics and Fatigue 

Westergaard [1939], Irwin [1957], and Williams [1957] were the first to write 
closed form solutions for the stress distribution near a flaw. Their solutions were 
limited to very specific geometries and loading conditions. Their results, in the form 
of a series solution, showed that the stress a distance r from a crack tip varied as r ~''' . 
It can be shown that, under linear elastic conditions, the first term of the series solution 
for the stress near a flaw in any body, under mode I, or opening, loading is given by: 

V 

where r and Q are polar coordinates as defined in Figure 3.1,/iy is a function of d that is 
dependent on the mode of loading, and K, is the mode I stress intensity factor. The 
sub- and super-scripts (/) denote mode I loading. Similarly, two other modes of 
loading can be defined as in-plane shear, mode II, and out-of-plane shear, mode III. 
The stress solutions for mode II and III loading are identical in form to Equation (3.1), 
but with all of the sub- and super-scripts I replaced with II or III. 
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Figure 3.1 : Coordinate system at a crack tip. 


A significant feature of Equation (3.1) is that as r goes to zero, or as one 
approaches the crack tip, this first term of the series solution approaches infinity. 
However, the higher order terms of the series will remain finite. For this reason, a 
large portion of LEFM focuses on this first term of the series expansion only. In 
reality, the stresses do not approach infinity at the crack tip. There is a zone around 
the tip where linear elastic conditions do not hold and plastic deformation takes place. 
This zone is called the plastic zone and results in blunting of the sharp crack tip. 
LEFM holds when the plastic zone is small in relation to the length scale of the crack. 

TTie SIF is a conveni^t way to describe the stress and displacement 
distributions near a flaw in Hnear elastic bodies. The SIF for any mode is a function of 
geometry, crack length, and loading. The general equation for a SIF is 

K = (3.2) 

f5 is a dimensionless factor that depends on geometry, 2a is the crack length, and o is 
the far field stress. It can be seen from Equation (3.2) that the units of K are 

stress * -y/iength . 

For a crack to propagate, the energy supplied to the system must be greater 
than or equal to the energy necessary for new surface formation. When supplying 
energy to the system, the energy can primarily go into plastic deformation or new 
surface formation. LEFM assumes that all of the energy supplied goes into forming 
new surfaces. As a result, LEFM predicts the material at a crack tip will fail when the 
mode I SIF, A), reaches a critical intensity called the fracture toughness, Kic- Fracture 
toughness is a material property and by definition is not dependent on geometry. 
Therefore, the criterion for fracture, or crack propagation, under LEFM, in mode I, is 
Kj > . Standard tests can be performed to measure values of fracture toughness 

[ASTM 1997]. The tests subject a standard specimen to pure mode I loading. The 
crack growth direction under pure mode I loading is self-similar. In other words, the 
crack tip in Figure 3.1 under only mode I loading will extend along the .v-axis. 

However, it is rare that a crack is subjected to pure mode I loading. More 
realistically, the loading will be a combination of all the modes. The mixed mode 
loading affects the fracture criterion and crack trajectory. For example. Mode II 


NAS A/CR— 2000-2 1 0062 


22 



loading will turn, or kink, the crack away from self-similar crack propagation. There 
are several proposed methods to predict the direction of crack growth under mode I 
and II loading. The most widely accepted methods are the maximum principal stress 
theory [Erdogan et al. 1963], the maximum energy release rate theory [Nuismer 1975], 
and the minimum strain energy density theory [Sih 1974]. Due to ease of 
implementation and demonstrated accuracy, the maximum principal stress theory will 
be used in this thesis. The method is based on two assumptions. First, the crack will 
propagate radially from the crack tip. The second is that the crack will propagate in a 
direction that is perpendicular to the maximum tangential stress. In other words, the 
crack will kink at an angle 6,„ where Gee is a maximum. For mode I and II loading, 
assuming plane strain conditions, Ggg is 

7 i=cos— AT, cos" —- — AT,/ sin 0] (3.3) 

^|27D■ 2l 2 2 \ 

The direction of crack growth can also be shown to correspond to the principal stress 
direction. Setting the partial derivative of Ggg with respect to 6 equal to zero, the angle 
0„, will be that which satisfies the equation 

/T/ sin0 4- AT,, (3cos0 -l)= 0 (3.4) 

From Equation (3.4), it is seen that if K,i equals zero, i.e. pure mode I loading, 
then the crack will propagate at an angle equal to zero. Figure 3.2 illustrates 
schematically the angle of crack trajectory, 0„„ with respect to the crack front 
coordinate system. 



Figure 3.2: Angle of crack trajectory with respect to crack tip. 


3.2.1 Fatigue 

Cracks have been known to grow when the mode I SIF is less than Kic- In 
these instances, the flaw has been subjected to cyclic loading. Cyclic loading can 
produce fatigue crack growth at loads significantly smaller than the fracture toughness 
of the material. Figure 3.3 illustrates how cyclic loading is characterized by the tensile 
load range, AS, and the load ratio, R. R is defined as the ratio of minimum stress, S„,„„ 
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to maximum stress, 5„,a.v» which, due to similitude, is equal to the ratio of minimum 
mode I SIF, K„^in, to maximum mode I SIF, K^ax- 


c 

min 




C 

^ max ^ max 


(3.5) 


Cyclic load histories can also be classified as proportional or non-proportional. When 
the ratio of Ku to Ki is constant during the loading cycle, the loading is proportional. 
Non-proportional is the case when this ratio varies with time. 


Stress or SIF I 



There are three regimes of fatigue crack growth as demonstrated in Figure 3.4. 
Regime I is related to crack initiation and small crack effects. As noted on the plot, 
there is a threshold value, A^m, below which fatigue crack growth will not occur. For 
AISI 9310 steel, values for A^,y, are reported to range from approximately 3.5 ksi*in° ^ 
~ 12 ksi*in°’^ [Binder et al. 1980], [Forman et al. 1984], [Propjietory source 1998]. 
As the stress ratio goes from positive to negative, the threshold value" increases. 

Regime II is commonly referred to as the Paris regime. The work of this thesis 
will only focus on crack growth in regime II. Crack initiation, small crack effects, and 
unstable crack growth (regime HI) will be ignored. A seminal development in 
predicting fatigue crack growth was from [Paris et al. 1961] and [Paris et al. 1963]. 
They discovered that a crack grows in fatigue at a rate that is a function of AA"/. They 
proposed that the nature of the curve in regime II could be described by: 

^ = C{AK,r (3.6) 

clN 

where N is the number of cycles, and C and n were proposed as material constants. 
Equation (3.6) is commonly referred to as the Paris model. When the crack growth 
rate in regime II is plotted on a log-log scale as a function of AAT, the slope of the 
curve is n. If the curve is extrapolated to the vertical axis, the intercept is C. 
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In regime III, the crack growth is unstable. A crack can grow in fatigue only 
when K, < K . As a result, regime III is bounded on the right by AAT/c- 



Figure 3.4: Typical shape of a fatigue crack growth rate plot. 

Crack growth in regime II creates striations on the fracture surface in certain 
materials under appropriate loading conditions. It has been shown that the spacing 
between striations is roughly equal to the macroscopic crack growth rate daJdN 
[Forsyth 1962]. In general, ductile alloys, e.g. aluminum alloys, form the most well 
developed striations. The material of interest in this thesis, AISI 9310 steel, is capable 
of forming striations [Bhattacharyya et gl. 1979] [Au et al. 1981] [McElvily et al. 
1996]. Au et al. successfully correlated fatigue crack growth rates to fatigue striations 
in AISI 9310 steel. 

Paris first proposed C as a material property. However, experimental research 
has found that C varies as a function of the stress ratio. The crack growth rate 
increases as the stress ratio increases. Fatigue crack growth data in regime II from 
tests conducted at different stress ratios, plots as shown in the left graph of Figure 3.5. 
The spread in the curves is explained by fatigue crack closure [Fiber 1971]. In 
general, it has been found that a crack will prematurely close prior to the tensile load 
^ing entirely removed. The level of stress at which this premature closing occurs is 
Sop (or, due to similitude, Kop). Incorporating fatigue crack closure phenomenon into 
Paris’ model should collapse the curves into a single line (the right graph of Figure 
3.5). This is accomplished by plotting on the abscissa ~^op^-’ 

rather than AK. This single curve is referred to as the “intrinsic” fatigue crack growth 
rate. More details of fatigue crack closure will be discussed in Chapter 4. 
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Figure 3.5: Schematic of fatigue crack growth rate data in Paris regime at different 
stress ratios collapsing into a single “intrinsic” curve. 


Using Paris’ model, the amount of crack growth per cycle for a given cracked 
object and load history can be predicted from the SIFs. In computational fracture 
mechanics, the FEM or BEM is used to calculate the SIFs. Several ways to calculate 
SIFs using numerical methods include the displacement correlation method [Chan et 
al. 1970], stiffness derivative [Parks 1974], /-integral [Rice 1968], and the universal 
crack closure integral [Singh et al. 1998]. The displacement correlation technique is 
used in this work because it relies only on displacement information on the boundary 
near the crack tip and because the method is computationally efficient. The numerical 
analyses of the spiral bevel pinion are conducted using the boundary element method, 
which solves for displacement information only on the boundaries. The displacement 
correlation method is computationally efficient since only a single numerical analysis 
is adequate to calculate the SIFs, unlike some of the other techniques that require two. 
Additionally, the mode I, II, and HI SIFs are all calculated by the same method. 

The displacement correlation method utilizes the fact that the displacements 
near a crack tip are proportional to the SIFs. Under pure mode I loading, the opening 
displacement, , is given by [Owen et al. 1983] 


where 


M.. = 


K = 


Kl 

4/u \ 2 k 
3-v 


(2k + 1) sin 


1 

-sm — 

I I 2 JJ 


1 + v 


for plane stress 


(3.7) 


?c = 3 - 4v for plane strain 

is the shear modulus of the material, v is Poisson’s ratio, and $ is the angle between 
the location of the displacement ah<h the normal to the crack tip. Equation (3 .7) can be 
rearranged to solve for K / = f (u ^.) . Along the crack front 6 = 180° . Knowing the 

material properties (E (elastic modulus) and v), and the crack opening displacement 
My, at a given distance r from the crack front, K/ can be calculated. 


K,= 



8(l-v^ 


(3.8) 
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Similarly, equations for K,i and Km can be written as a function of «>, the displacement 
due to in plane shear, and m^, the displacement due to out of plane shear, respectively. 
It is important to note that as r approaches zero, the accuracy of the SIFs will decrease 
when using the displacement correlation method if the crack front elements are not 
capable of representing the singularities at the crack tip. 

Crack growth rates are calculated from the SIF information and experimentally 
determined fatigue crack growth model parameters. The SIF information is also used 
to calculate the angle of propagation, e.g. Equation (3.4). 

3.2.2 Example: Two dimensional, mode I dominant fatigue crack growth 
simulation with static, proportional loading 

The purpose of this example is to demonstrate how fatigue crack predictions 
can be performed on a simple two dimensional model. The model assumptions are: 

1 . The location of applied load is not changing. This will be referred to as 
static loading. 

2. The loading is proportional. 

3. The crack growth can primarily be attributed to mode I opening. In 
other words, K / » K ^ . This will be referred to as mode I dominant. 

4. Crack closure effects will be ignored. 

5. LEFM holds. 

The method to predict crack trajectories in two dimensions is incremental. A 
series of finite element analyses are run which incrementally increase the crack length 
by a significant amount in relation to the model’s geometry. For a given increase in 
crack length, the number of cycles to achieve that amount of growth can be calculated. 
For a given propagation step there are A, load cycles associated with it. 

The amount of crack growth for one cycle is calculated as a function of the 
maximum stress in the load cycle. Because it is assumed the loading is proportional, it 
is straightforward to calculate the direction the crack will grow during the cycle using 
the maximum principal stress theory. However, there is no proposed method to 
calculate the final amount and direction of crack growth during one load cycle if the 
ratio KulKi varies during the cycle, i.e. non-proportional loading. 
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Create 

geometry model 



Figure 3.6: Flow chart of process to predict fatigue crack trajectory. 

As outlined in Figure 3.6, the process begins with a geometry model. The 
geometry model is then discretized into a finite element mesh. Figure 3.7 shows the 
finite element mesh for an arbitrary geometry model that will be used for 
demonstrative purposes. This particular initial mesh consists entirely of quadratic 
eight-noded elements. 

Model attributes must be defined next. The material properties are specified 
within the finite element program as a Young’s Modulus of 29,000 ksi and Poisson’s 
ratio of 0.25. The thickness of the model is taken to be 1 inch. Boundary conditions 
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must also be defined. A cyclic loading history like that shown in Figure 3.3 is 
assumed. The minimum applied traction is assumed to be zero, and the maximum 
applied traction is = 100 ksi . A tensile traction is applied normal to the top edge. 
All of the nodes along the bottom edge are restrained in the vertical direction, and the 
far right comer node is restrained in the horizontal direction. If desired, at this stage 
the finite element solver could be mn to calculate displacements, strains, and stresses 
in the uncracked geometry. 



Figure 3.7; Two dimensional finite element model. 

Next, a crack is introduced into the geometry model. With the change in 
geometry, the model must be remeshed. However, the damage to the mesh model is 
localized, and, therefore, only a small region around the crack must be remeshed. The 
mesh around the crack tip is a rosette of eight triangular, six-noded quarter point 
elements. Figure 3.8a. The remaining area is meshed with quadratic six-noded 
elements. Figure 3.8b shows the initial edge crack and locally remeshed region. The 
boundary conditions, material properties, and loads were defined earlier and do not 
need to be redefined. At this point, displacements, strains, and stresses are solved for 
in the cracked geometry. 

A method, such as the displacement correlation technique, is used to compute 
the stress intensity factors at the crack tip based on the relative displacements of the 
crack faces. Once the SIFs are calculated, Paris’ model (Equation (3.6)) can be used 
to calculate the amount of crack growth for one load cycle, daldN. A method, e.g. 
maximum principal stress (Equation (3.4)), is used to determine the direction of crack 
growth from the calculated SIFs. In most cases, the amount of crack growth for one 
load cycle will be on the order of 10^ - 10^ inches. Since this is significantly smaller 
than the geometry features of the gear, it would be inefficient to update the geometry 
model for every load cycle. Consequently, a number of load cycles is assumed, e.g. 
Nj = 2,000 cycles . Finally, the crack in the geometry model is extended by an amount 
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daj = (da / dN)i *A^y, at an angle 6„, with respect to the self-similar crack trajectory 
(Figure 3.2). Again, the model must be remeshed locally, and the process is repeated. 




a) Initial crack and quarter point element b) Finite element mesh after adding initial 
rosette crack 

Figure 3.8: Initial edge crack in model. 



Figure 3.9: Predicted crack trajectory for model in Figure 3.8. 


Figure 3.9 is a picture of the predicted crack trajectory for the finite element 
model in Figure 3.8. The crack has been incrementally advanced from the initial 
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length and orientation through five propagation steps. For the assumed material 
properties and loading in this example, the calculated number of load cycles to grow 
the crack from the initial length in Figure 3.8 to that in Figure 3.9 is 4,900 cycles. 

3.2.3 Example: Three dimensional, mode I dominant fatigue crack growth 
simulation with static, proportional loading 

The assumptions of the two dimensional example will apply to this three 
dimensional example. In three dimensions, the procedure to predict fatigue crack 
trajectories is very similar to that in two dimensions. As in two dimensions, the 
geometry model must be defined, the mesh created, and the model attributes assigned. 
The main complexity with three dimensional crack growth simulations is that there is 
not a single crack tip, but rather a three dimensional crack front. For a given three 
dimensional crack, there is no longer a single value for the SIF in each mode, but 
rather a SIF distribution along the crack front for each mode. In addition, the crack 
length might also vary along the crack front. 

In this thesis, all of the three dimensional models are boundary element 
models. In the boundary element method, the primary variables are load and 
displacement. Strains and stresses are secondary variables. The BEM is based on an 
integral equation formulation. An advantage of the method is that the number of 
unknowns in the equations is proportional to the surface discretization. This is in 
contrast to the FEM where the number of unknowns is proportional to the volume 
discretization. In computational fracture mechanics when predicting crack trajectories 
and remeshing are necessary, an advantage of the BEM is that only the surfaces near 
the crack need to be remeshed, as opposed to the entire volume which must be 
remeshed when using the FEM. Volume meshing with cracks can be rather difficult; 
whereas, surface meshes are straightforward with and without cracks. 


crack face 
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There are no closed form solutions to calculate SIF distributions along the 
crack front for arbitrary three dimensional cracks. As a result, a conventional 
approach to calculate the SIF distribution is to discretize the front into a series of two 
dimensional crack tips. For example, the finite plate model presented in Section 3.2.2, 
in reality, has a finite width. Therefore, the crack must have a finite width. The crack 
front shape might be that shown in Figure 3.10. In this example, the crack width is 
equal to the plate thickness. 



Discretized three dimensional crack front Two dimensional crack tip 

Figure 3.1 1 : Discrete crack front points treated as two dimensional problems. 


Next, the crack front is discretized, as shown by the lines intersecting the crack 
front in Figure 3.11. Once the crack front is discretized, each point is treated as a two 
dimensional problem. The two dimensional methods to calculate SBFs are applied at 
each discrete point. The discrete point is propagated by an amount and at an angle 
uniquely defined by the SIFs associated with that point. Once each discrete crack 
front point is propagated individually, a least squares curve fit is performed through 
the new discrete crack front points. Figure 3.12. _ 

A potential difference in the three dimensional approach, as opposed to the two 
dimensional method, is that singular crack front elements might not be used along the 
crack front. Since the BEM is implemented in this thesis, the volume of the three 
dimensional model is not meshed; only the surfaces are meshed. Therefore, elements 
that represent the crack tip singularity are not available along the crack front. The 
main drawback of this is that some SIF accuracy along the crack front is sacrificed. 
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New discrete crack 
front points 


Figure 3.12: Least squares curve fit through new discrete crack front points. 

3.3 Fracture Mechanics Software 

A suite of fracture mechanics software developed by the Cornell Fracture 
Group is used in this thesis [FRANC3D 1999a, 1999b]. The codes were developed to 
handle the complexities of three dimensional crack trajectory predictions. OSM is 
used to define a three dimensional solid geometry model of an object. The program is 
based on defining the surfaces of the model explicitly in Cartesian space. The 
boundary of a solid is generated by adjacent surfaces, or faces. Each face of the 
boundary element model has a three dimensional local coordinate system associated 
with it. In order to define a closed solid, all of the local face normals must point away 
from the interior of the solid. The local coordinate system might also be of 
significance when defining boundary conditions. 

The geometry model is then read into FRANC3D. With FRANC3D, a user 
can create a finite element or boundary element mesh based on the geometry model. 
Displacement or force/traction boundary conditions must be defined for all the faces 
of the solid. The conditions must be specified in all three Cartesian directions with 
respect to either the local or the global coordinate system. Material properties are also 
assigned to regions of the model using FRANC3D. 

Cracks are added to the solid by explicitly defining the vertices, edges, and 
faces that model the cracks. A crack has two distinct faces that must be meshed 
identically. 

As mentioned in Section 3.2.3, a crack front must be discretized prior to 
calculating SIFs and to propagating the crack. Within FRANC3D, there are three 
options to discretize the crack front. The discrete points can be defined by the mesh 
nodes, the midpoints of the elements sides along the crack, or at a user defined number 
of equally spaced points along the crack front. The built in feature in FRANC3D to 
calculate SIFs uses the displacement correlation technique. The most accurate results 
are obtained when a row of four sided elements is used along the crack front. This 
will give a set of equally spaced points behind the crack front where the SIFs can be 
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evaluated. Additionally, to improve the performance of the crack front elements, the 
ratio of the elements’ width to length should be close to one [FRANC3D 1999c], 

When a crack is propagated, the geometry model changes. However, the 
geometry changes only near the crack. Therefore, only the mesh model near the new 
crack is damaged and requires remeshing. The remainder of the geometry and mesh 
model is left unchanged. This is a distinct advantage of FRANC3D. 

The program BES is used to solve for the displacements and stresses using the 
boundary element technique. FRANC3D is used as a post-processor to view the 
deformed shape, stress contours, and extract nodal information. 

FRANC3D uses the same functional form to interpolate the geometry and field 
variable variations over an element. The form is given by the associated element type. 
In all of the models, only isoparametric three- and four-noded elements are used. 
Quadratic elements are available; however, based on the work in [FRANC3D 1999c], 
the gain in accuracy does not justify the significant increase in computational time. 

3.4 Chapter Summary 

This chapter covered theories of LEFM and fatigue pertinent to modeling crack 
growth numerically. Of primary importance is how crack growth rates and trajectory 
angles are calculated from SIFs. The maximum principal stress theory will be used to 
calculate trajectories under mixed mode loading. In addition, the displacement 
correlation method was introduced as a technique to evaluate SIFs. Two dimensional 
and three dimensional examples demonstrated how the theories are applied in 
numerical simulations. Some features of the software programs FRANC3D, OSM, 
and BES that will be used in the simulations were covered. The background provided 
In Chapter 2 and this chapter will be utilized in the work of Chapters 4, 5, and 6. The 
studies in those chapters cover issues related to predicting three dimensional fatigue 
crack trajectories in a spiral bevel gear. 
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CHAPTER FOUR: 

FATIGUE CRACK GROWTH RATES 

4.1 Introduction 

The goal of this chapter is to determine how highly negative stress ratios affect 
the fatigue crack growth rate in a common gear steel, AISI 9310. This is of interest in 
the context of gears because the magnitude of compressive stresses in a gear’s tooth 
root is a function of the rim thickness. If fatigue crack growth rates are highly 
sensitive to compression, then crack growth rates may warrant more attention in 
designing gears. On the other hand, if the compressive stresses do not alter the fatigue 
crack growth rate predictions greatly, than the loading cycle for a gear tooth can be 
simplified by ignoring the compressive portion of the cycle. 

In Section 4.2, the concept of fatigue crack closure is discussed. This section 
shows that crack closure provides a convenient framework within which to understand 
the factors that control fatigue crack growth. A material-independent method is 
presented for obtaining fatigue crack growth rate data that does not vary with stress 
ratio, R. The crack closure approach is extended beyond aluminum alloys, considered 
by Elber [1971] and Newman [1981] and discussed in Section 4.2, to steels. Next, 
Section 4.3 applies the concepts to AISI 9310 data to obtain an intrinsic fatigue crack 
growth model. Section 4.4 demonstrates that in the range of negative R, the effective 
stress range, and likewise the crack growth rate, is not highly sensitive to the 
magnitude of R. 

4.2 Fatigue Crack Closure Concept 

Due to the cyclical loading on a gear’s tooth, fatigue crack propagation might 
occur. The load range, A5, or stress intensity factor range, A^, along with the load 
ratio, R, characterizes cyclic loading. Recall, R is defined as the ratio of minimum 
stress, Smin, to maximum stress, Smax, which, due to similitude, is equal to the ratio of 
minimum mode I SIF, to maximum mode I SIF, K„,a,x (Equation 3.5). 

Lewicki et al. [1997b] found that spur gear teeth can have /?-values as low as 
-3.0. They also found that the magnitude of /? in spur gears is a function of the gear 
geometry. As the rim thickness decreases, R becomes more negative due to the 
increased bending of the gear rim. 

A general interpretation of the crack closure approach is that damage only 
occurs during the portion of the load cycle when the crack faces are not in contact. 
TTie majority of the literature’s discussion of crack closure covers its effect on crack 
growth rates. Since gears have such high load frequencies, crack growth rates are 
commonly of secondary interest in the context of gears. The time from detectable 
flaw to failure is usually insignificant. However, if the crack growth rate is highly 
sensitive to the magnitude of the compressive portion of the load cycle, then crack 
growth rates may warrant more attention. On the other hand, if, for negative values of 
R, the crack growth rate is relatively insensitive to the magnitude of R, then the effect 
of geometry on R need not be the primary concern in gear design. This demonstration 
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is significant in the context of the overall goal of this thesis, which is to study aspects 
of gear geometry that affect damage tolerance. 

It is assumed initially in this chapter that the stresses induced in a gear tooth 
under positive (tensile) and negative (compressive) parts of the load cycle are 
“proportional.” In other words, the shape of the stress intensity factor distribution 
along the crack front is the same under both tensile and compressive loading. In two 
dimensional analyses, this is not a concern because the crack only consists of a tip, 
where the deformation can be tensile only or compressive only, not a combination of 
the two. In three dimensions, however, the distribution of the loading (deformation) 
along the crack front might be different in the compressive and tensile load cases. In 
the end, whether the positive and negative parts of the load cycle are proportional is 
not of major concern. As will be shown in the remaining sections, damage occurs 
only during the tensile portion of the load cycle. 

Elber [1971] observed that during unloading a crack actually closes prior to the 
applied load being entirely removed. This phenomenon has been called fatigue crack 
closure. Fatigue crack closure also explains why, for a given AK, fatigue tests show 
the crack growth rate increasing as R increases. Figure 4.1 shows typical fatigue crack 
growth rate data as a function of SIF range [Kurihara et al. 1986]. Kurihara et al. 
conducted fatigue tests with 500 MPa class C-Mn steel, which is used in pressure 
vessels. The tests covered a wide range of stress ratios from -5.0 to 0.8. Figure 4.1 
was obtained by selecting two data points off Kurihara et a/.’s plots for each value of 
R. The horizontal scatter in the curves is a result of the different /(-values. Note that 
as R increases, the curves shift to the left, producing an increase in fatigue crack 
growth rate for a given AK. 
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Figure 4.1: Fatigue crack growth rate data for pressure vessel steel at various /(-values 

(data taken from [Kurihara et al. 1986]). 
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Crack closure can be attributed to a number of factors. During the opening 
portion of a load cycle, the material at the crack tip plastically deforms. As the cycles 
repeat, a wake of plastic deformation remains as the crack propagates through the 
body. The plastic deformation wake results in a mismatch between the crack faces. 
Although not considered here, crack closure can also occur due to differences in the 
surface roughness of the crack faces, due to mixed mode loading, or oxidation of the 
crack surfaces. 

Elber modified Paris’ model to account for crack closure. The modification 
allows crack propagation to occur only while the crack tip is open. He introduced Sop 
as the stress level where the crack first opens during the tensile part of the load cycle. 
His equation for the crack propagation rate is: 

^ = C(AK^r =C(UAKy (4.1) 


where U, the effective stress range ratio, is defined as 


S - S 

max op 




\-R 


(4.2) 


Figure 4.2 illustrates the relationships among various K values. Sop {Kop) is difficult to 
measure experimentally. In addition, the value varies with loading conditions. As a 
result, Elber developed an empirical relationship between U and R. From this 
relationship. Sop (Kop) could be backed out. 



Figure 4.2: Constant for different stress ratios. 

When daldN is plotted as a function of AKeff, the scattered curves (due to 
different /^-values) collapse into a single, “intrinsic” crack growth rate curve. In 
crack-closure-based fatigue models, daldN is a function of AKcg. This implies that 
crack growth occurs only while the crack tip is open. If Kmax is kept constant between 
various tests with different 7?-values, then Ar„„„ must change. If it can be shown that 
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remains nearly constant for various negative i?-values, then the portion of the 
load cycle when < K^p does not contribute to crack growth. Therefore, all 

negative 7?-value cases could be treated in the same manner. The sensitivity of ^Keff 
will be investigated in Section 4.4. Figure 4.2 illustrates how AK^ff could remain 
constant as Kmm decreases. 

Fiber performed a series of experimental fatigue tests with sheets of 2024-T3 
aluminum alloy. The stress ratio range was -0.1</?<0.7. From the tests, he 
developed the empirical relationship 

U(R) = 0.5 + OAR when - 0. 1 < ^ < 0.7 (4.3) 

Fiber’s U(R) relationship is valid only for 2024-T3 aluminum alloy over the 
range of /?-values for which he had experimental data. His work inspired many to 
develop empirical relationships between U and R for a variety of materials and ranges 
of R. Schijve [1988] summarizes several of these empirical relationships for different 
alloys and ranges of R. However, it is expensive to develop this relationship 
empirically every time one wants to model crack closure in a new material. This led 
to attempts to numerically model crack closure [Newman 1976, 1981], [Fleck et al. 
1988], [McClung et al. 1989], and [Blom et al. 1985]. Through the thirty plus years 
of research related to crack closure, it has been found that the amount of crack closure 
is dependent on many variables. Specimen size, specimen geometry, crack length, 
applied stress state, and prior loading conditions all affect the magnitude of Sop- 
Newman’s work attempts to incorporate all of these factors. 

Newman developed and applied a hybrid analytical/numerical crack closure 
model that simulated plane strain and plane stress conditions. He successfully 
matched crack growth rates under constant-amplitude loading from his analytical 
model to experimental data. The material he focused on initially was 2219-T851 
aluminum alloy. The model has since been applied to a variety of metals. Newman’s 
model is the most comprehensive and has been successfully validated against 
experiments. As a result, his model will be utilized in this thesis. 

All variables in Fquation (4.2) are defined immediately from the loading conditions 
with the exception of Sop- To find an expression for SopISmax, Newman [1984] fit 
equations to his numerical results for 2043-T3 aluminum alloy over a large range of R- 
values and load levels. He worked in terms of applied loads, but due to similitude, S 
in his expressions can be replaced with K, giving; 

^^ = A, + A,R + A,R^^+A,R^ for R>0 (4.4a) 

= A, + A.R for - 1 < 7? < 0 (4.4b) 


when K^,p > . The coefficients Ao - Aj are: 


= (0.825 - 0.341C + O.OSk' ) 


COS 


V 2(To I 


(4.5) 


NASA/CR— 2000-2 10062 


38 


A, =(0.415;c-0.071>c-)- 


A2 — 1 -^0 ^3 

A3 = 2 Aq + a, — 1 

vis a constraint factor taking on a lower bound of 1 for plane stress conditions and an 
upper bound value of 3 to simulate plane strain conditions. The flow stress Ob, is 
avLage between the uniaxial yield stress and the uniaxial ultimate tensile strength 

the material. ^ ^ i? 

Because Newman’s model for is a function of material constants (ao), K, 

and K, it is applicable for any fatigue crack where LEFM holds and 
conditions and material properties are known. Figure 4.3 is an examp e 0 0 

curves in Figure 4.1 collapse into an intrinsic curve when crack closure is taken 1 
account. Equations (4.4) and (4.5) are used to calculate AKop- U is calculated using 

Equation (4.2).^ 



The crack tip condition in the fatigue test specimen Kunhara ct al. used, a thin 
plate with a center crack, is best described by plane stress. Therefore a value of 
V = 1 was selected for the preliminary graphs. K was then mcreased, and the amount 


5 Note that Newman claims Equation (4.4b) is valid for negative /?-values greater • 

H^ver. Kurihara er a/.’s data extends to -5.0. Equation (^4b) was used 

-3.0, and -2.0. Figure 4.3 illustrates, at least for this case, the equation can also hold for these low R 


values. 
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of correlation between the curves was visually inspected. As k increased, the curves 
became more scattered, validating the choice of jc = 1 . 

The equation of a line in Figure 4.3 is given by: 

In Figures 4.1 and 4.3, the slope for a given curve (/?-value) is uniquely defined by the 
data points. According to the crack growth models, all of the curves should have the 
same slope. Ideally, this would be the case for the plots in Figure 4.3. The relatively 
small scatter in the magnitude of the slopes at different /?-values is attributed to scatter 
in the experimental results. 

Figure 4.3 includes the intrinsic curve predicted by the intrinsic ^ = 0 data. 
This curve falls roughly in the middle of the predicted curves. To give an idea of the 
scatter in the curves, the figure also includes lines corresponding to one half and two 
times the crack growth rate for /? = 0. All of the predicted intrinsic curves fall into 

this envelope. As a result, it is concluded that the AK^ff equations produce good 
correlation. 

These results with 500 MPa pressure vessel steel demonstrate that an intrinsic 
fatigue crack growth rate curve can be obtained using Newman’s material-independent 
model to account for crack closure. It is also shown that a possibility exists to extend 
the model beyond the range of /?>-!. Consequently, in Section 4.3 the model will 
be applied to AISI 9310 steel to determine how negative /?- values influence crack 
propagation rates. 

4.3 Application of Newman’s Model to AISI 9310 Steel 

An open literature search for fatigue crack growth rate data for AISI 9310 steel 
at various /?-values revealed little published information. A report by Au et al. [1981] 
contains the most information. Au et al. performed tests in different environments at 
various /?-values and frequencies for carburized and noncarburized AISI 9310 steel. 
Because they were investigating the correlation between fatigue striations and crack 
growth rates, only two tests were performed on noncarburized steel in the same 
environment and at the same load frequency but at different I?-values. The load levels 
used in the tests were not reported. When their measured fatigue crack growth rates at 
R = 0.05 and 0.5 are plotted against dX, there is very little scatter in the curves. This 
suggests that the crack growth rate is not sensitive to R or that the applied loads were 
high enough such that Since the objective of this study is to correlate 

fatigue crack growth rate data at different /?-values, including the negative R regime, 
Au et o/.’s data is inadequate. 

Additional fatigue test data for AISI 9310 was provided by a helicopter 
manufacturer on the condition that the data’s source not be identified. Data points are 
extracted from the fatigue crack growth rate curves obtained from tests in two 
different environments. Figure 4.4 shows growth rates for AISI 9310 steel in room 
temperature air for /?=-!, 0.05, and 0.5. The curves in Figure 4.5 are obtained by 
extracting data points from fatigue crack growth rate tests in 250° oil for /? = -1, 0.01, 
and 0.5. Table 4. 1 summarizes the slopes and intercepts for the various curves. 
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Figure 4.4: Fatigue crack growth rate data for AISI 9310 steel in room temperature air. 
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Figure 4.5: Fatigue crack growth rate data for AISI 9310 steel in 250"" oil. 
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Table 4.1: Slope and intercept of curves in Figures 4.4 and 4.5. 


Test n 

C 

f ( in! cycle}! ( ksi* irP'^)^] 

R = A (Air) 

3.3 

6.4e-12 

R = 0.05 (Air) 

3.5 

7.3e-ll 

R = 0.5 (Air) 

3.9 

5.3e-ll 

R = -l (Oil) 

3.2 

l.le-11 

I? = 0.01 (Oil) 

3.2 

9.9e-ll 

R = 0.5 (Oil) 

3.8 

7.9e-ll 


For a given R, the n values are similar between the two environments. The 
effect of the environment can be see in the variations of C. C is consistently larger in 
the heated oil environment. A larger C will result in faster growth rates. However, the 
environment effect will not be considered in this investigation. 

Similar to the pressure vessel steel analyses, da/dN versus AKeff plots are 
generated using Equations (4.2), (4.4), and (4.5). A value of tc = 1 best describes the 
condition at the crack tip in the test specimen. Figures 4.6 and 4.7 illustrate the 
various curves collapsing into an intrinsic fatigue crack growth rate curve. 


l.OOE-03 , 



AATgj [ksi*in®'’] 

Figure 4.6: Intrinsic fatigue crack growth rate for AISI 9310 in room temperature air, 

k: = 1. 

Figures 4.4 through 4.7 demonstrate that Newman’s crack closure model 
accounts for the scatter in fatigue crack growth rates at different stress ratios in AISI 
9310 steel. Table 4.2 contains the slopes and vertical intercepts from the lines in the 
figures. In addition, a linear least squares curve is fit through the data in Figures 4.6 
and 4.7. The slope and vertical intercept from each curve fit are also included in the 
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table. These curve fit values will be used in the crack growth rate models for the 
numerical analyses. Crack closure concepts will now be extended to investigate the 
sensitivity of fatigue crack growth rates to low /? -values. 
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Figure 4.7: Intrinsic fatigue crack growth rate for AISI 9310 in 250° oil; X" = 1 . 


Table 4.2: Intrinsic and non-intrinsic fatigue crack growth rate parameters. 


Test 

C 

** [( inlcycle)l(ksi*irf'^f'] 

c* 

[ ( in! cycle)! ( ksi*uf')"] 

R = -l (Air) 

3.3 

6.4e-12 

6.30e-10 


3.5 

7.3e-ll 

8.80e-10 


3.9 

5.3e-ll 

1.98e-10 

R = -l (Oil) 

3.2 

l.le-11 

8.52e-10 

/? = 0.01 (Oil) 

3.2 

9.9e-ll 

1.09e-9 

I? = 0.5 (Oil) 

3.8 

7.9e-ll 

2.87e-10 

Curve Fit Air 

3.6 

NA^ 

4.26e-10 

Curve Fit Oil 

3.4 

NA* 

7.44e-10 


6 


Not Applicable 
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4.4 Sensitivity of Growth Rate to Low R 

Table 4.3 contains results from calculations of AKeff Ht different /?-values using 
Equations (4.1), (4.2), (4.4), and (4.5). Constant values for K„,ax, K and SmaJoa are 
assumed. 


Table 4.3: Calculations to find AAT^j^over a range of /?-values for a constant K,nax\ SIF 

units are ksi*in”^. 


R 

K 

Stnax^ ^0 

f^max 


AK 

U 

AKeff 

0.705 

1.000 

0.100 

10.000 

10.000 

7.050 

7.000 

2.950 

0.825 

2.434 

0.700 

1.000 

0.100 

3.000 

0.822 

2.467 

0.505 

1.000 

0.100 

10.000 

5.050 

4.950 

0.716 

3.542 

0.500 

1.000 

0.100 

10.000 

5.000 

5.000 

0.713 

3.565 

0.255 

1.000 

0.100 

10.000 

2.550 

7.450 

0.589 

4.388 

0.250 

1.000 

0.100 

10.000 

2.500 

7.500 

0.587 

4.399 

0.005 

1.000 

0.100 

10.000 

0.050 

9.950 

0.474 

4.714 

0.000 

1.000 

0.100 

10.000 

0.000 

10.000 

0.472 

4.716 

-0.495 

1.000_ 

0.100 

10.000 

-4.950 

14.950 

0.327 

4.886 

-0.500 

1.000 J 

0.100 

10.000 

-5.5o 

15.000 

0.326 

4.888 

-0.955 

1.000 

0.100 

10.000 

-9.950 

19.950 

0.254 

5.058 

”-1.000 

1.000 

0.100 

10.000 

-10.000 

20.000 

0.253 

5.060 

-1.995 

1.000 

0.100 

10.000 

-19.950 

29.950 

0.180 

5.402 

-2.000 

1.000 

0.100 

10.000 

-20.000 

30.000 

0.180 

5.404 

-2.995 

1.000 

0.100 

10.000 

-29.950 

39.950 

0.144 

5.746 

-3.000 

1.000 

0.100 

10.000 

-30.000 

40.000 

0.144 

5.748 


The crack growth rate is calculated in Table 4.4 based on the effective SIF data 
in Table 4.3. C and n are assumed to be 7.44e-10 (in/cycle)/(ksi*in”'^)" and 3.4, 
respectively. These values are taken from the curve fit to the intrinsic growth rate data 
for the AISI 9310 steel tests conducted in heated oil. daldN as a function of R is 
plotted in Figure 4.8. 

The curve in Figure 4.8 shows that the crack growth rate in the negative R 
regime is less sensitive to variations in R compared to the positive R regime. Between 
R equal to zero and -3.00, the crack growth rate varies by a factor of 1 .96. In a fatigue 
context, a difference of this order of magnitude is acceptable. As a result, one can 
conclude that when modeling fatigue crack growth, AKe^, or likewise Kop or daldN, 
does not change significantly for R<Q . Therefore, the magnitude of R is not a useful 
parameter to characterize damage evolution in gears. In the context of designing gear 
geometry to be damage tolerant, a primary concern need not be how aspects of gear 
geometry affect R. 

It has been shown that when crack closure is taken into account there is not a 
significant change in the crack growth rates for negative /?-values. This result will be 
utilized in the numerical analyses discussed in Chapters 5 and 7. The load ratio will 
be taken as = 0 under the assumption that, if /? < 0 , the general results and 
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conclusions would still be valid. This is a simplification to the loading cycle and 
method. 

Table 4.4: Crack growth rate calculations for a wide range of /?-values taking into 
account crack closure effects. The percent change in daldN is due to AR = 0.005 . 


R 

[ksi*in°'^] 

da/dN 

(in! cycle x 
10'^) 

% Change 
daldN 

0.705 

2.434 

0.153 

— 

0.700 

2.467 

0.160 

4.446 

0.505 

3.542 

0.549 

— 

0.500 

3.565 

0.561 

2.135 

0.255 

4.388 

1.136 

— 

0.250 

4.399 

1.146 

0.873 

0.005 

4.714 

1.1449 

— 

0.000 

4.716 

1.451 

0.130 

-0.495 

4.886 

1.637 

— 

-0.500 

4.888 

1.639 

0.120 

-0.955 

5.058 

1.841 

— 

-1.000 

5.060 

1.844 

0.116 

-1.995 

5.402 

2.303 

— 

-2.000 

5.404 

2.306 

0.108 

-2.995 

5.746 

2.841 

— 

-3.000 

5.748 

2.844 

0.102 



Figure 4.8: daldN as a function of AKgff. 
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4.5 Chapter Summary 

Highlights of this chapter can be summarized as follows: 

• The crack closure concept and Newman’s model were presented. It was shown 
that the model predicts that fatigue damage occurs only during the portion of the 
load cycle when the crack faces are not in contact. 

• Newman’s crack closure model was applied to empirical data for crack growth 
rates of a pressure vessel steel. It was shown that crack closure explains well the 
apparent dependence of crack growth rates on R. In fact, the material has an 
intrinsic crack growth rate. /? is a parameter that determines during what portion 
of the load range the crack faces are not in contact. This range is called the 
effective stress intensity factor range. 

• Newman’s model was applied to AISI 9310 steel, a typical steel used for gears. 
There was much less crack growth data available for this steel as compared to the 
pressure vessel steel. Nevertheless, it was shown that the crack closure model 
works for this small data set. 

• It was demonstrated that, in the regime of negative f?-values, the model predicts 
that the crack growth rates as a function of the effective stress intensity factors are 
only a weak function of the magnitude of /?. 

• The observation made in this chapter that crack growth rates are not highly 
sensitive to R in the negative 7?-regime will be used in Chapter 5 when modeling 
the load history. 
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CHAPTER FIVE: 

PREDICTING FATIGUE CRACK GROWTH 
TRAJECTORIES IN THREE DIMENSIONS UNDER 
MOVING, NON-PROPORTIONAL LOADS 

5.1 Introduction 

Chapter 5 covers numerical modeling issues related to predicting fatigue crack 
growth trajectories in three dimensions in a spiral bevel pinion. The goal of this 
chapter is to model crack growth under realistic operating conditions. As covered in 
Section 1 .2, most previous work in the area of predicting crack trajectories in gears 
assumed one fixed load location. The location was usually the HPSTC. However, in 
operation, spiral bevel gears are subjected to a load moving in three dimensions. The 
fixed location loading, therefore, could lead to incorrect three dimensional trajectories. 

A boundary element model of the OH-58 spiral bevel pinion is presented in 
Section 5.2. The tooth coordinates and a dimensioned drawing of the pinion were 
provided by NASA/GRC, along with the coordinates for discrete elliptical contact 
areas along a spiral bevel gear tooth. OSM/FRANC3D is used to create the model 
from these data. Studies are conducted to determine the smallest model that still 
achieves accurate SIF results. Once this model is defined, initial analyses for the 
discrete load cases are conducted. The SIF history for an initial crack subjected to the 
moving load is presented in Section 5.3. 

Section 5.4 develops a method to predict three dimensional fatigue crack 
growth trajectories under a moving load. The method increments a set of discrete 
points along the crack front for each step in the load cycle to find the total amount of 
extension and final angle of growth after fifteen load steps (1 load cycle). The 
propagation path for each point is then approximated with a straight line. A number of 
cycles are specified, and the crack front is advanced an amount equal to the crack 
extension for one cycle times the number of assumed cycles, and at the angle 
calculated for one cycle. Next, a curve is fit through the new crack tip locations to 
define the new crack front. The FRANC3D geometry model is updated, and the 
process is repeated. 

Finally, in Section 5.5 the proposed moving load crack propagation method is 
implemented to predict fatigue crack growth trajectories in the OH-58 spiral bevel 
pinion. 

5.2 BEM Model 

A boundary element model of the OH-58 spiral bevel pinion was built with 
OSM/FRANC3D. The Cartesian coordinates for a tooth surface, tooth profile, and 
fillet curve were provided by NASA/GRC. The data were generated automatically 
from a program that models the gear cutting process along with the gear kinematics. 
All points on the generated tooth surface are points of tangency to the cutter surface 
during the manufacturing process [Litvin 1991]. A primary motivation for developing 
the tooth geometry program was to generate data for a three dimensional finite 
element analysis. This program’s output was adapted to develop a boundary element 
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model for this thesis. The remainder of the pinion solid model was built from a 
drawing of the pinion. The basic shape of the shafts and gear rim were modeled. 
Some subtle details of the pinion drawing were ignored in cases where the geometry 
would complicate the geometry model and have negligible effects on the computed 
SIFs. The surfaces of the solid model were meshed using three- and four-noded 
elements. Figure 5.1 contains three views of a typical boundary element model (recall 
that the meshes shown in the figures are surface meshes). The volume of the gear is 
not meshed. The conical shape of the gear rim and the cylindrical shape of the shafts 
are seen best in Figure 5.1b. As seen in Figures 5.1a and 5.1c, three of the nineteen 
teeth of the pinion are modeled explicitly. Section 5.2.2 discusses studies to verify the 
accuracy of the three teeth model. 



a) Overall view of full model 



Gear rim 


b) Section A- A from (a): profile of shaft 
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c) Close up view of teeth 

Figure 5.1: Typical boundary element model of OH-58 pinion. 

In operation, the input torque is applied at the end of the pinion’s long shaft. 
The small shaft sits on roller bearings. When the torque is applied, the gear rotates 
and the teeth of the pinion successively contact the gear’s teeth. When contact occurs, 
load is transferred across the teeth. The boundary conditions shown in Figure 5.1a 
model these operating conditions. This model will be referred to as the full model. 
The face patches at the end of the long shaft are fixed in all directions. The 
displacements on the surfaces of the smaller shaft are restrained in the local normal 
direction. Though not explicitly shown in Figure 5.1, contact areas are modeled as 
distinct face patches on the middle tooth. Traction normal to the patch is defined 
which equals the load that is transferred across the contacting teeth for a given input 
torque and rotation angle. More detail on how these contact patches are defined is 
given in the next section. 

5.2.1 Loading Simplifications 

The meshing of the gear and pinion is a continuous process. The magnitude of 
force between the gear teeth varies during the meshing as adjacent teeth come into and 
out of contact. Figure 2.8 is a schematic of the continuous process that has been 
discretized into fifteen load steps. 

In order to perform numerical crack propagation studies of the pinion, the 
continuous contact between the teeth is discretized into fifteen contact patches, or load 
steps: four double tooth contact patches, seven single tooth contact patches, followed 
by four more double tooth contact patches. Each load step is a unique face patch in 
the boundary element model. The load steps will be referred to as numbers one 
through fifteen, corresponding to the patches from the gear root to the top land, 
respectively. This is consistent with the progression of contact area along a pinion 
tooth from the root toward the top. One progression through the fifteen load steps is 
one load cycle on the tooth. One rotation of the gear results in one load cycle on each 
tooth. 
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The location and size of the fifteen discrete contact patches were provided by 
NASA/GRC. The data were determined numerically by a procedure similar to that 
described by Litvin ef al. [1991]. The mean point of contact between the gears is 
taken as the center of the ellipse. Hertzian contact theory along with the applied 
torque level is used to determine the width of the ellipse. The patches were calculated 
for operating conditions of 300 horsepower, 6060 rotations per minutes, and 3 1 20 in- 
lb torque. These conditions are approximating the 100% design load condition, which 
is defined as 3099 in-lb torque. 



Figure 5.2: Contact ellipses defined at the geometry level in the numerical models. 
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In the BEM model, the shape of a contact ellipse is approximated by straight 
lines connecting the axes' end points. The straight line approximation is valid because 
Saint Venant’s principle holds; as long as the total applied forces and resulting 
moments are kept constant, the elliptical shape of the traction can be approximated by 
a patch with straight sides without altering the stress distribution along the crack front. 
Frictional forces are neglected, and, consequently, the traction is constant over the 
patch. Each patch has a unique magnitude of traction. 

The four figures in Figure 5.2 demonstrate how the traction patches are built 
into the model geometry. The purpose of the models is to calculate SIFs from all 
fifteen static load cases. The combination of all fifteen SIF distributions represents 
one load cycle on the tooth. Figure 5.2 shows how a single BEM model can 
incorporate multiple load cases. Not all of the contact ellipses can be modeled in one 
BEM model because there is overlap between the ellipses. The multiple load case 
feature minimizes the computational time. For example, the boundary element model 
for load cases one, five, eight, and thirteen is virtually identical. The only difference 
between them is the boundary conditions. Hence, with the multiple load case feature, 
the two most computationally expensive steps of the boundary element solver, setting 
up the boundary integral equations and factoring the stiffness matrix, occur only once. 
The different boundary conditions are then applied individually, and the corresponding 
equations are solved for the unknown displacements and tractions for each load case. 

5.2.2 Influence of Model Size on SIF Accuracy 

The fewer the number of elements, or unknowns, in a boundary element 
model, the less computationally intensive the model is. Minimizing the number of 
elements can primarily be accomplished by 1) using a coarser mesh with larger 
elements or 2) by modeling less of the geometry of the solid. A disadvantage of the 
first option is the accuracy of the solution is sacrificed. The elements used in all of the 
studies in this thesis are linear. Therefore, only linear variations in displacement 
across an element can be represented. Likewise, the geometry is approximated by a 
series of linear segments. Because the geometry of the pinion is complex with 
significant amounts of curvature, larger elements do not represent the geometry 
adequately. As a result, this option is disregarded, and the second option, simplifying 
the model, is considered. 

Simplifying the model also has drawbacks. The smaller the portion of the 
pinion modeled, the less accurate the representation of the boundary conditions. Three 
simplified models are investigated. The first simplification. Figure 5.3, is to ignore the 
long shaft in the full model. The new faces that are created when the long shaft is 
disregarded are restrained in all directions. Secondly, the smaller shaft is removed. 
Figure 5.4. The boundary conditions on the heel end are the same as simplification 
one, and the new faces on the toe end are set to traction free. The final simplification 
is to cut the rim of the pinion in half. Figure 5.5. The boundary conditions for this 
model are the same as the second simplification, with the addition of roller boundary 
conditions (displacement in the direction of the local normal set to zero) applied to the 
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new faces. The boundary conditions for each model are chosen because they most 
closely match those of the full model (Figure 5.1a). 



Figure 5.3; Simplified model I: ignore long shaft. 


In each of the simplified models, the flexibility of the pinion changes. When 
an identical crack is introduced into all of the models, the SIFs might vary from model 
to model. To determine whether a simplified model is valid, the SIFs from the 
simplified models are compared to the full model’s SIFs for identical cracks. It is 
assumed that the full model most accurately represents the operating conditions and 
loading paths. 
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Figure 5.5: Simplified model 3: ignore both shafts and half of gear rim. 

A semi-elliptical crack is introduced in the root of the middle tooth in each 
model. The crack is 0.64 inches long and 0.14 inches deep. A simplified load is 
applied to the middle tooth over the middle third of the tooth length and tooth height. 
The shape of the traction patch is rectangular, and the traction across the patch is 
constant. The shape and location of the traction patch is different from those 
described in Section 5.2.1. However, the difference is not important because the intent 
is to analyze differences in SIFs between models after changing one variable and 
keeping all the rest of the model parameters constant. To achieve consistency between 
all the models, the mesh in the region of the crack and load patch is identical. 

The SIFs increase on average by 7%, 8%, and 11% with respect to the full 
model’s SIFs for simplification one, two, and three, respectively. In a fatigue growth 
rate context, changes of this magnitude are significant. Recall that the crack growth 
rate is proportional to Ki raised to a power (Equation (3.6)). For AISI 9310, the 
magnitude of the exponent is approximately 3.4. Consequently, seemingly small 
changes in the SIFs have dramatic effects on the crack growth rate predictions. It is 
concluded from this study that the full model should be used for all trajectory 
predictions. 

To verify that only explicitly modeling three teeth yields accurate results, a 
nine teeth model is analyzed. If the SIFs between the three teeth and nine teeth 
models are similar, then it can be concluded that not all of the nineteen teeth of the 
pinion need to be modeled. 

An edge crack is introduced in the three and nine teeth models, in the middle 
of the tooth length, in the root of the concave side of the middle tooth. The crack 
shape is semi-elliptical, and is 0.125 inches long and 0.05 inches deep. An effort is 
made to keep the meshes between the two models identical. 

The difference in SIF distribution under load steps 1, 5, and 8 is investigated. 
As shown in Figure 5.6, the percent difference in Ki between the two models for all 
three load cases is below 5%. The.absplute magnitude of Ku for both models and all 
load cases is significantly smaller than Ki. Consequently, a small variation in Ku 
appears as a large percent difference between the models. Instead of percent 
differences. Figure 5.7 shows the absolute AT// values for all the loads and models. It is 
evident from the figures that the three teeth and nine teeth models yield similar results. 
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leading to the conclusion that the three teeth model is sufficient for the trajectory 
prediction analyses. 



-6 4 


-8 - 
-10 

Figure 5.6: Percent difference in Ki between three teeth and nine teeth models for load 
cases one, five, and eight. Crack front position one corresponds to the heel end of the 

crack front. 



Figure 5.7: Kn distribution for three teeth and nine teeth models for load cases one, 

five, and eight. 

5.3 Initial SIF History Under Moving Load 

To simulate the moving load during one load cycle on a pinion’s tooth, fifteen 
static BEM analyses are performed. Each analysis represents ope of the fifteen 
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discrete time steps as the contact area moves up the pinion tooth, as discussed in 
Section 5.2.1. Recall these contact ellipses are defined for a full design load input 
torque of 3120 in-lb. The full pinion boundary element model is used. Figure 5.1. A 
semi-elliptical edge crack is introduced into the root of the middle tooth on the 
concave side. The crack is located approximately in the middle of the tooth length. 
The dimensions are 0.125 inches long by 0.050 inches deep. The crack is oriented 
approximately normal to the surface. 

Each load step produces a unique SIF distribution along the crack front. The 
SIF distribution changes between load steps because the load position and magnitude 
varies from step to step. Figure 5.8 shows the mode I SIF distribution for the first 
eleven load steps, the initial four double tooth contact load steps followed by the seven 
single tooth contact load steps. The second stage of double tooth contact, load steps 
twelve through fifteen, are omitted from the figure to simplify it. Modes II and HI 
SIFs are plotted similarly in Appendix A. 

The four bottom curves in Figure 5.8 are the SIFs under double tooth contact. 
The remaining seven curves are the SIFs under the single tooth contact load steps. 
The bottom most of the seven curves corresponds to load five. The topmost curve is 
the result from load eleven, the last single tooth contact step. The total applied force 
for each single tooth contact load step is roughly equivalent. However, as the load 
step number increases, the SIF curves shift up. This is explained by the fact that the 
locations of the contact patches are progressing up the pinion tooth. The change in 
location creates a greater moment arm. As a result, the displacements, and likewise 
SIFs, in the tooth root will thus increase. 



21 26 31 36 41 

Crack front position 
(Orientation: heel to toe) 


-^Load 1 
-■-Load 2 
-A- Load 3 
-H- Load 4 
— Load 5 
Load 6 
— Load 7 
— — Load 8 

Load 9 

-•-Load 10 
-B-Load II 


Figure 5.8: Mode I SIF distribution for load steps one through eleven. 


NAS A/CR— 2000-2 1 0062 


55 



Another approach to examine the data is to plot the SIF history for each point 
along the crack front. Figure 5.9 shows the SIF history for point 29 in Figure 5.8 
(roughly the midpoint of the crack front). The magnitude of AT/, AT//, and K,n is plotted 
as a function of time, or load step. The figure also includes Kop, which was calculated 
using Newman’s crack closure equations described by Equations (4.1), (4.2), (4.4), 
and (4.5). 


10000 1 



When the individual points in Figure 5.9 are connected with straight lines, the 
plots represent the loading cycle on the tooth. The minimum load has been taken to be 
zero. In actuality, the minimum load in the tooth root might be compressive. When a 
tooth is loaded, compressive stresses could result in the root of the convex side. 
Depending on the magnitude of these stresses, they may extend into the concave root 
of the adjacent tooth. However, Chapter 4 demonstrated that the crack growth rates do 
not vary significantly for negative /?-values when crack closure is taken into account. 
Therefore, the load cycle is modeled as 7? equals zero, i.e. AT/„„„ = A^ = K =0. 

The difference in the single tooth and double tooth contact loads is evident in 
the mode I SiFs. The plateaus in the curve correspond to the two contact stages. Kj is 
significantly larger during single tooth contact (load steps 5-11) compared to the 
double tooth contact stages (load steps 1-4 and 12-15). TTie magnitudes of AT/ are 
significantly greater than AT///. As a result, it will be assumed that mode in does not 
contribute to the crack growth. 

Based on gear theory, the curves in Figure 5.9 should be continuous and 
smooth. The continuous curves would most likely show a large increase in slope as 
the loading progresses from double tooth contact to single tooth contact. One can 
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imagine that as the number of discrete load steps increases, the curves in Figure 5.9 
will become smoother. However, due to transmission error and noise, it is known that 
the curves in reality are neither continuous nor smooth. Therefore, it is assumed that 
the fifteen load steps are sufficient to approximate the true loading conditions. 

The moving load on the pinion’s tooth is non-proportional; the ratio of Ku to Ki 
changes during the load cycle, Figure 5.10. Consequently, the method to propagate a 
three dimensional crack described in Section 3.2.3 can not be used. That method 
assumed proportional loading, which results in a constant kink angle for the load 
cycle. Since the ratio is changing in the spiral-bevel gear tooth, the predicted angle of 
propagation during the load cycle changes. A method to determine how a crack would 
grow under this type of loading is required and is proposed in Section 5.4. 



The Ku to Ki ratio also indicates which loading mode is driving the crack 
growth. Mode I dominant fatigue crack growth is associated with smaller ratios. Qian 
et al. [1996] studied mixed mode I and II crack growth in four point bend specimens. 
They selected the test specimen geometries from FEM analyses that considered 
various crack lengths and orientations to achieve different Ku to Ki ratios. From the 
analyses, they selected five different geometries with KulKi values of 0, 0.262, 0.701, 
1.812, and 16.725. The ratios covered crack growth mechanisms of pure mode I, 
mode I dominant, balanced mode I and II effects, mode II dominant, and highly mode 
n dominant, respectively. Using these ratios as guidelines, the gear situation can be 
characterized as balanced mode I and n effects during the earlier stages of the cycle to 
mode I dominant crack growth during the later stages of the cycle. However, it will be 
assumed that the fatigue crack growth is driven by mode I. 
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5.4 Method for Three Dimensional Fatigue Crack Growth Predictions 
Under Non-Proportional Loading 

As shown in Section 5.3, a crack in a spiral bevel pinion tooth is subjected to 
non-proportional loading. As a result, conventional methods to predict crack growth 
trajectories in three dimensions are not adequate. A literature review of non- 
proportional fatigue crack growth revealed only a few methods that were applicable to 
the gear model; Section 5.4.1 is a summary of relevant work. A method to predict 
three dimensional fatigue crack trajectories under non-proportional loads is proposed 
in Section 5.4.2. Section 5.4.3 summarizes the approximations of the method. 

5.4.1 Literature Review 

In the literature related to non-proportional fatigue crack growth, the majority 
of the work is experimental. The limited amount of numerical work is related to 
predicting crack growth rates and fatigue life. The numerical work is also largely 
confined to two dimensional analyses. Schijve [1996] gives an overview of methods 
and research related to predicting fatigue life and crack growth. There is no mention 
of predicting crack trajectories in non-proportional loading scenarios. Crack 
trajectories are of primary importance in the context of gears because the trajectory 
determines whether the failure mode will be catastrophic. The number of cycles to 
failure is of secondary importance because the high loading frequency on a gear’s 
tooth results in very short times from crack initiation to tooth or rim failure. 

Bold et al. [1992] is the most extensive report covering fatigue crack growth in 
steels under mixed mode I and II loading. They give experimental results from non- 
proportional mode I and II tests, and compare the maximum tangential stress theory 
(pure mode I) and maximum shear stress (pure mode II) theory for predicting kink 
angles to experimental results. However, their work contains no theoretical 
predictions for mixed mode. Bower et al. [1994] considered brittle fracture under a 
moving contact load. They incrementally advanced the load and evaluated the SIFs at 
each stage of contact. If the SIFs met their fracture criterion, then the crack was 
propagated based on the mode I and II SIFs for that load position using the maximum 
principal stress criterion. Their approach incorporates non-proportional loading in an 
incremental manner; however, the work is limited to brittle fracture, does not include 
fatigue, and does not include three dimensional effects. 

Hourlier et al.’s [1985] focus was to determine which of three theories 
predicted trajectories closest to experimental data for non-proportional loading. They 
worked in terms of A/, which is the mode I stress intensity factor for a small advance of 
the crack at an angle 9. The three theories investigated were 1) direction in which A/ is 
a maximum, 2) direction where AA/ is a maximum, and 3) direction of maximum 
fatigue growth rate da/dN. The rate is calculated assuming a mode I dominant growth 
mechanism and is a function of ki„,a.x(6) and Akj(6). Their work found that 1) was the 
most inaccurate method. In general, the maximum daldN method was found to best 
match experimental results. Hourlier et al.'s work is not practical for the purposes of 
this thesis primarily because it requires one to express the moving load and A/ in closed 
forms as functions of time and 6 in order to find the angle corresponding to the 
maximum daldN. 
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Three dimensional finite element analyses have been performed to simulate the 
wheel position over a railroad track containing a crack [Olzak et al. 1993]. The rail 
model is analyzed for consecutive stages of wheel position and the SIFs are calculated 
for each stage. However, Olzak et al. did not propagate the crack. Their primary goal 
was to determine what happens to the crack displacement and contact shape when the 
load is directly over the crack. In the case of the spiral bevel gear, the load will most 
likely never be directly over the crack and their findings are not applicable. 

The most significant work was done by Panasyuk et al. [1995]. They 
numerically modeled and propagated a two dimensional edge crack under a moving 
contact load. The maximum principal stress theory was used, and growth rates were 
calculated by Paris’ model. The translation and location of the contact are expressed 
as functions of A, the distance from the load to the crack. To calculate the kink angle, 
first the values of A that correspond to an extremum of A' = F[A'/(A), AT//(A), 0(A)] are 
found. Next 6, Ki, and ^// at these A are calculated, from which the growth rate is 
calculated. Finally, it is assumed that the crack propagates for N cycles at that growth 
rate and angle, and the crack in the numerical model is updated and the process is 
repeated. Panasyuk et al. assumed that their geometry was an elastic half plane, and, 
therefore, they could set up closed form equations and solve analytically for A'/, A'//, 
and 6. Once again, their method can not be directly applied to gear model because 
neither the traction nor the geometry can be expressed in closed form. The method 
also does not directly take into account the non-proportional loading and assumes a 
constant kink angle for the entire load cycle. However, their method is extended and a 
similar incremental approach is developed in the next section. 

5.4.2 Proposed Method 

Compared to a two dimensional static problem, the problem at hand is 
continuous in time and in a third space dimension. Methods have been presented in 
previous chapters and sections to discretize both of these dimensions. With the 
discretizations, two dimensional crack propagation theories can be applied. In 
summary, the proposed method discretizes the continuous loading in time into a series 
of elliptical contact patches, of load increments. Two dimensional fatigue crack 
propagation theories are then used to propagate incrementally a series of discrete 
points from the three dimensional crack front. The remainder of this section outlines a 
proposed method to predict fatigue crack trajectories in three dimensions taking into 
account time varying SIFs. 

Method 

1 . Discretize tooth contact path into 15 load steps (Section 5.2.1). 

2. Calculate by the displacement correlation method, using a feature built in to 
FRANC3D/BES, the mode I, II, and El SIFs (A^^j, Kuq), Kn/ij)), where i is a 
discrete point along the crack front (i=l-num _points) and j is the load case 
(j = 0-M). In general, the nodes of the first row of mesh nodes behind the crack 
front are taken as the discrete points. Figure 5.9 is a typical plot of the SIFs for a 
single point i along the crack front for the entire loading cycle. 
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3 . 


The goal of this step is to calculate for a given point / the amount of crack 
extension, during a load step from j-\ to j. 6'(j.\j) is the angle for the 

extension during a load step from j-l to j. It is assumed that the crack grows 
incrementally during a load cycle. In addition, propagation at point i only takes 
place when the change in mode I SIF between load steps is positive, 

(Kj' j -K/ij-i)) >0, and only when K/j is greater than the opening SIF at that 
point, Kop'. This implies that growth will only take place during the loading 
portion of the cycle. 

To calculate most accurately the total amount of crack growth over one cycle, 
crack closure is taken into account. The amount of extension during one load 
cycle predicted by a modified Paris’ model, adjusted to incorporate crack closure, 
is 

da‘=C(AK,J)'' (5.1) 

where AKeJ = K/„,ax - Kop = U is given by Equation (4.2). Kop is found 

using Equations (4.4) and (4.5). Figure 5.9 shows that the loading is characterized 
by ^ = 0 . 

In order to calculate U, S„ax, the far field stress in a Griffith crack problem, is 
required. Figure 5.11 shows the Griffith crack geometry [Griffith 1921]. The 
gear’s geometry is obviously different from a Griffith crack problem. Therefore, 
an equivalent Smax must be calculated for the gear. First, K/„,ax is found in step 2. 
Smax is then found by solving Equation (3.2) for Smax'' 



Lastly, it is assumed that, at a given point, the amount of extension between 
load steps is proportional to the ratio of the change in mode I SIF to the effective 
SIF. The amount of crack growth for each load increment is given by: 

K,\jy-K,'u-u 


da\j-i.j) — 


AK 


-da' 


(5.3) 




The angle of crack growth associate with each load increment is found from the 
maximum principal stress theory using the current load steps SIFs as: 


e\ 


7 - 1 . j) 


= 2 tan' 


1 K,\j, ^ 1 


^K' 


It 


f 

I (y> 


K 


It <>> 


+ 8 


(5.4) 
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Figure 5.11: Griffith crack problem: Straight, through thickness crack in an infinite 
plate subjected to uniform tensile stresses [Griffith 1921]. 


4. 


Repeat step 3 for every load step of the cycle to get the final coordinates and angle 
for the trajectory during one load cycle. The final crack trajectory is approximated 
by a straight line from the initial crack tip location to the final crack growth 
location. Based on simple geometry, the final length, daj, and final angle. Of, 
after one load cycle are calculated in the following manner: 

da,' +{h'y- (5.5) 


where 



= tan 





) 


/' = ^da'(j-uj)COs(p'(j-ij)) 
;=i 

h = ^ Ja'o-I.;) sin((9'o-i.y) j 

;=i 


(5.6) 

(5.7) 

(5.8) 


Figure 5.12 illustrates this step schematically, assuming the load cycle has 
been discretized into four steps, i.e. M = 4. Note that the arc length, which is the 
sum of the day.\j), is equal to the amount of growth predicted by the crack- 
closure-modified Paris model. Equation (5.1). The arc length is given 
mathematically by 

da' ~'^da'(j-\,j) (5.9) 

j=0 

5. Repeat steps 3 and 4 for every point along the crack front. 
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6. Determine the number of cycles, N, necessary to achieve a significant amount of 
crack growth in relation to the pinion’s geometry. This step is necessary because 
the amount of crack growth over one load cycle is too small to update the 
geometry model. Therefore, it is necessary to assume that a series of load cycles 
has occurred prior to changing the geometry. Because each point along the crack 
front has a unique growth rate associated with it, the crack front will not grow 
uniformly. Each point will grow by an amount daf,„ai- 

(5-10) 

In general, N is chosen such that dap„a/ > 0.01 inches. 

7. Update the FRANC3D geometry model with the new crack that has grown by an 
amount of dafi„J- To accomplish this, a least squares curve fit is performed 
through the new discrete crack front points. A single polynomial curve may be fit 
through all of the points, or the points may be divided into a user-defined number 
of sets and individual polynomial curves are fit to each set. After the crack 
geometry is updated, it is necessary to locally remesh the model prior running the 
BEM solver. Once again, all load steps are analyzed with the new crack. 

8. Repeat process beginning at step 2. 



Figure 5.12: Schematic of crack extension for one point along the crack front after one 

load cycle. 
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5.4.3 Approximations of Method 

The accuracy of the proposed method is limited by the accuracy of the BEM 
results. In addition, several approximations and simplifications in the loading 
conditions and crack growth rate models could affect the accuracy. The shape of the 
elliptical contact patches is estimated by straight lines. When the applied traction is 
far enough away from the crack front such that the deformation due to local contact 
stresses does not contribute to the crack tip field, St. Venant’s principle holds. In this 
case, the shape of the traction area is inconsequential and the straight line 
approximation is legitimate. However, if the crack trajectory is significantly close to a 
contact patch, then this assumption is no longer valid. As a result, the accuracy of the 
SIFs could be comprised. Furthermore, the size and location of each of the fifteen 
contact areas are kept constant throughout the crack propagation. In reality, the 
change in flexibility of the tooth could change the contact areas and in turn affect the 
SIF distribution in the later stages of propagation. This aspect will be investigated in 
Chapter 7. However, if the same contact areas are used during the entire crack 
propagation simulation, the fatigue life predictions will most likely be conservative. 
This is because the SIFs continuously increase as the crack advances when the loading 
scenario is kept constant and the crack length is increasing. 

It is assumed that the crack front conditions are characterized by plane strain 
along the entire front. This is consistent with the maximum principal stress theory 
equations implemented to calculate the kink angle. However, shallow qracks and 
portions of edge cracks near the free surface are usually characterized by plane stress, 
not plane strain. Nevertheless, crack growth rates will be larger in plane strain 
conditions than plane stress. This assumption errs on the conservative side and will 
predict shorter fatigue lives. The size of the plastic zone in the gear is investigated 
more in Section 7.3.2. 

An additional approximation is introduced when incrementally calculating the 
amount of crack growth. The method assumes that crack growth only occurs during 
the tensile portion of the load cycle. However, it can not be experimentally shown that 
this is true. In fact, it is generally accepted that crack growth occurs during the 
opening and closing portion of the load cycle [Laird 1967]. The tensile portion creates 
plastic deformation at the crack tip and causes it to blunt. During unloading, the 
plastic deformation creates a wedging action at the crack tip that acts to advance the 
crack. 

The reasoning behind Equation (5.3) governing the crack growth rate is that, if 
the loading were to become proportional, the amount of crack extension predicted 
during one loading cycle would be equal to that predicted by the crack-closure- 
modified Paris’ model. An additional assumption is that the method assumes mode I 
dominant fatigue crack growth. If the ratio of Kn to AT/, and likewise the ratio of Km to 
AT/, were to become large enough, the mode II (or mode III) loading could contribute to 
the crack growth. In this case, should be a function of Kj, Ku, and/or Km- The 
models proposed for crack growth rates do not incorporate the mode II and HI effects. 

For each propagation step, a value for N is chosen that is large enough to 
increase a majority of the crack front a significant distance in relation to the model’s 
geometry. It is assumed during the N loading cycles that the variations in the 
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displacements and stresses are negligible. In reality, the crack is growing slowly 
during the cycles and therefore the displacements and stresses at the crack front will 
vary as the crack advances. If the value of N is too large, then the model will not pick 
up subtle changes in displacements. Consequently, the predicted crack front shape 
and trajectory will be inaccurate. 

5.5 Simulation Results 

Three dimensional fatigue crack trajectory simulations are performed using the 
method described for non-proportional loading in Section 5.4.2. The full pinion 
boundary element model is used and an initial crack is introduced into the model as 
described in Section 5.3. In order to validate the predictions, the dimensions and 
shape of the crack are similar to a starter notch that was used in a tested pinion. The 
experimental results and data are reported in Chapter 6. 

Thirteen crack propagation steps are carried out, which model the crack at 
thirteen distinct stages of growth. In each propagation step, the SIFs are evaluated 
using the displacement correlation method at the mesh nodes in the first row of 
elements along the crack front. To reduce the computational time for each 
propagation step, the second stage of double tooth contact, load steps twelve through 
fifteen, are ignored. The method for propagation adopted here only allows crack 
growth during loading. The final four load steps represent the unloading portion of the 
load cycle and, therefore, do not contribute to crack growth in the simulations. 

The Paris growth rate model modified to incorporate crack closure is used in 
this study. The model parameters are held constant during the propagation steps. 
Values for C and n are taken from a curve fit to AISI 9310 steel crack growth rate data 
from 250° oil (Section 4.4). K is set equal to three (plane strain). To calculate Sma.x. 
which is a variable in the Kop calculation, /3 is set to one. The material properties listed 
in Table 2.2 are used. 

A least squares curve fit to the predicted discrete crack front points is used to 
determine a smooth crack front curve. The approach is to fit a polynomial curve of 
second or third order to groups of points. To allow unsymmetric crack front shapes, 
the crack front points are divided into one, two, or three groups. A curve is then fit 
independently through each group. 

The tooth contact locations and magnitudes defined for the 100% design load 
are used throughout the thirteen propagation steps. The loading simplifications for 
these data were presented in Section 5.2.1. 

Table 5.1 contains the crack geometry and growth rate data for each of the 
propagation steps. N is rounded to the nearest 100 cycles. Figures 5.13 and 5.14 show 
the initial and final crack trajectories along the tooth surface and the depth of the crack 
into the gear rim. At the end of propagation step thirteen, the crack has propagated 
neither entirely through the tooth width nor through the rim. The analyses are stopped 
at this step because the toe end of the crack has reached the top land. It was decided 
that continuing the analyses would lead to no additional insights because, as reported 
in Chapter 7, the prediction differed from the experimental results. Also, one can 
imagine that if the simulations were continued, since the trajectory has turned on both 
the heel and toe end, this middle portion of the tooth will break away from the gear. 
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Propagation 

Step 

N 

[cycles] 

Crack Front 
Length [in] 

Depth^ 
[in X 10'"] 

Area 

[in- X 10'-] 

0 

N.A. 

0.200 

5.00 

0.579 

1 

15,000 

0.233 

5.53 

0.743 

2 

20,200 

0.258 

5.58 

0.837 

3 

38,200 

0.237 

5.17 

0.933 

4 

56,200 

0.453 

5.88 

1.69 

5 

66,500 

0.506 

6.32 

2.19 

6 

76,800 

0.595 

6.85 

2.77 

7 

94,900 

0.774 

8.10 

4.27 

8 

121,000 

0.940 

9.50 

5.76 

9 

147,200 

1.11 

11.4 

7.77 

10 

187,400 

1.15 

12.9 

10.1 

11 

227,000 

1.19 

14.9 

12.5 

12 

274,000 

1.29 

16.9 

15.7 

13 

311,000 

1.42 

18.8 

18.6 

^Toial = 

311,000 






b) Cross section of tooth at midpoint of crack 
Figure 5.13: Initial crack; N = 0 cycles. 


^ The approximate location along the tooth length of the initial crack’s midpoint is used to measure the 
crack depth. 
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b) Cross section of tooth at midpoint of 
initial crack 

Figure 5.14: Crack prediction after thirteen propagation steps; N = 31 1,000 cycles. 

In the simulation, the majority of the crack front has progressed almost 
halfway through the tooth width. A rough estimate of the total life of the gear is made 
by estimating the tooth width as 0.227 inches at the toe end and using the average 
growth rate from step thirteen. From these data, approximately 200,000 cycles are 
necessary for the crack front to progress through the remaining ligament. Therefore, 
the number of cycles to failure beginning from the initial notch is estimated to be 
511,000. This number could be non-conservative since a constant K is used to 
calculate the remaining life. On the other hand, the fatigue life prediction does not 
take into account the cycles leading up to crack initiation. It assumes that the crack 
begins propagating immediately after the introduction of the notch. It is most likely 
that in a real gear a number of cycles are attributed to initiation of the crack 
propagation from the notch. 

Figure 5.15 is a plot of the calculated fatigue life. The location along the tooth 
length of the initial crack’s midpoint is used to measure the crack depth, which is 
plotted against total number of cycles. The crack face area as a function of total 
number of cycles is also plotted for comparison since the crack depth varies along the 
crack length. 


NAS A/CR— 2000-2 1 0062 


66 
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N [cycles] 

Figure 5.15: Predicted crack depth and crack area versus number of cycles. 

The accuracy of the crack trajectory and fatigue life predictions will be 
evaluated in Chapter 7 by comparing the simulation results to experimental data. The 
experimental data is obtained from a spiral bevel pinion test that was carried out by 
NASA/GRC. 

5.6 Chapter Summary 

In this chapter, a boundary element model of a spiral bevel pinion was 
presented. Different size models were investigated to determine the smallest model 
which achieved accurate SIF results. An initial crack was introduced into the model 
and the SIF history along the crack front under the moving load on the pinion tooth 
was found. It was determined that the loading on the tooth was non-proportional. As 
a result, a method was developed to propagate the three dimensional crack under the 
non-proportional load. A crack trajectory prediction in the OH-58 spiral bevel pinion 
was performed using this method. In Chapter 7, the simulation results will be 
compared to experiment results; Chapter 6 presents the experimental data. In addition, 
in Chapter 7 parameters of the crack propagation method will be investigated to 
determine the sensitivity of the growth rates and trajectories to variations in the 
method and crack growth rate model variables. 
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CHAPTER SIX: 
EXPERIMENTAL RESULTS 


6.1 Introduction 

Experimental results can be used to evaluate the accuracy of the predictions 
made in Chapter 5. Recordings of crack length and depth as a function of total 
number of cycles are necessary to confirm the fatigue life prediction. Knowing the 
crack front shape during propagation would also assist in verifying the predictions. 

There is a limited amount of useful experimental data for the OH-58 pinion to 
validate the numerical results. Test data from an OH-58 spiral bevel pinion are 
provided by NASA/GRC. The crack growth observations made during the test are 
limited. As a result, the fracture surfaces of the tested pinion are observed under a 
scanning electron microscope (SEM). The test data from NASA/GRC and the 
fractography results are summarized in Section 6.2 and 6.3, respectively. A crack 
growth scenario during the test is formulated from the SEM observations. 


6.2 Test Results 

A pinion that was tested by NASA/GRC in their spiral bevel gear test fixture 
under a separate research project is used for comparison/validation. Table 6.1 
contains the loading data from the test. A notch was electro-discharge machined 
(EDM) into the root of a tooth’s concave side. The gear was run for six million cycles 
beginning at 1550 in-lb torque and progressing up to 4649 in-lb torque at the end of 
the six million cycles. The test was not stopped until the completion of the six million 
cycles; at which time, there was no observable crack growth. As a result, eight more 
notches of varying sizes were fabricated into individual tooth roots of the pinion. The 
pinion then ran continuously for an additional 4.9 million cycles at increasing levels of 
torque detailed in Table 6.1. At the completion of the 4.9 million cycles, five teeth 
had fractured from the pinion. 


Table 6.1 : Pinion test data. 


Time 

[cycles] 

Speed 

[rpm] 

Torque 

[in-lb] 

# of EDM 
Notches 

1 million 

6060 

1550 

1 

1 million 

6060 

2324 

1 

1 million 

6060 

3099 

1 

1 million 

6060 

3874 

1 

2 million 

6060 

4649 

1 

1 million 

4848 

2479 

9 

1 million 

6060 

3099 

9 

1 million 

6060 

3874 

9 

1 .9 million 

6060 

4649 

9 


During the latter 4.9 million cycles, the test was never stopped to observe the 
crack growth. Therefore, the sequence of events for the tooth fractures and the exact 
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number of cycles to failure are unknown. Some of the crack growth may have 
occurred during the first million cycles at 2479 in-lb torque, or there may have been 
no crack growth until the last 1.9 million cycles at 4649 in-lb torque. In addition, 
since no observations of crack growth during the 4.9 million cycles were made, the 
predicted fatigue life curve (Figure 5.15) can not be validated. The only quantitative 
information is an upper bound of 4.9 million cycles on the total number of cycles to 
failure. 
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Qualitative information can be taken from the pinion test to validate the 
predictions. For example, all five failures were teeth fractures. There was no 
evidence of a rim type of failure. A macroscopic view of the crack trajectory can be 
determined from the tested pinion. Figure 6. 1 shows photographs of a typical tooth 
failure in the tested pinion. The rectangular EDM notch on the concave side of the 
tooth near the root is observable in all of the pictures. The crack trajectory, which is 
assumed to initiate from the EDM notch, is deeper into the rim in the middle of the 
tooth length than on the toe end. In addition, a ridge is observed where the crack 
trajectory turned toward the root of the convex side from the initial path into the rim. 
All five of the fractured teeth had remaining portions of the tooth left intact at the heel 
end. Furthermore, four of the notches were small enough such that they do not appear 
to have grown during the 4.9 million load cycles. 

Crack growth rate data and crack front shape information during propagation 
could both be used to validate the numerical analyses. Since this data was 
unavailable, the tooth fractures were observed with a SEM. The SEM observations 
are summarized in the next section. 

6.3 Fractography 

6.3.1 Overview 

Three teeth that had fractured from the tested OH-58 pinion were examined 
using a SEM. Tooth #1 1 contained a “short” notch (~4 mm long x ~lmm deep), and 
was the initial flaw from which the numerical simulations were based. Tooth #5 
contained a longer and deeper starter notch (~18 mm x -1.7 mm deep). Pictures of 
these teeth’s fracture surfaces are contained in this section. The starter notch of the 
third tooth examined was similar to tooth #11. No distinguishing features were found 
on the tooth that were not seen on the other teeth. 

There were several goals of the examinations. The primary objectives were to 
determine how much of the crack growth occurred due to fatigue, where the 
mechanism of failure might have changed, and, if it did change, to what type of 
fracture mechanism. Additionally, information on the rate of crack growth and the 
crack front shape during propagation was desired. All of this information could verify 
how well the numerical simulations predicted the actual crack evolution. 

The sketch in Figure 6.2 shows a typical view of how the crack propagated 
through a pinion tooth’s cross section. The crack initiated from the EDM notch on the 
side of the tooth with the applied load. This side is the loaded or concave side. The 
opposite side of the tooth, where the crack ends, is the load free or convex side. The 
region where the crack deviates from the original path, into the rim, towards the fillet 
on the convex side will be referred to as the ridge. 

Figures 6.3 and 6.4 are of tooth #5 and #11, respectively. The fractured teeth 
are lying on their side with the tooth surface at the bottom of the figure and the ridge 
at the top. As illustrated in Figure 6.2, the crack growth direction was from the bottom 
of the photograph toward the ridge at the top of the photograph. The EDM notch is 
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the rectangular region along a portion of the tooth length near the bottom right comer 
and is labeled A in the figures. 


Concave side Convex side 



— Fatigue: rubbed 

Fatigue: partially rubbed 

— ■ Ductile rupture 


Figure 6.2: Sketch of crack propagation through typical pinion tooth’s cross section. 



Figure 6.3: Low magnirication (4.6x) view of footh"#5^s fracture surface on concave 

side. 

The left side of both Figures 6.3 and 6.4 is the toe end of the tooth length. This 
is true for all of the concave side pictures in this chapter. On the other hand, the right 
side of the photographs from the convex side is the toe end of the tooth length. All of 


NAS A/CR— 2000-2 1 0062 


72 




the photographs of the fracture surfaces designate the tooth number, the level of 
magnification, and a length scale. 



Figure 6.4: Low magnification (5.3x) view of tooth #1 1 ’s fracture surface on concave 

side. 


6.3.2 Results 

The morphologies of teeth #5 and #11 are very similar. Lines emanating 
radially from the starter notch are visible at low magnifications (4.6x in Figure 6.3 and 
5.3x in Figure 6.4). These lines are indicative of fatigue crack growth. 

On tooth #11, some fatigue striations can be seen near the EDM notch. Region 
A in Figure 6.5a is the EDM notch. Figure 6.5b is a magnified view of the striated 
region.^ As expected, the striations are roughly parallel to the edge of the EDM notch. 
The crack growth direction was perpendicular to the striations and, as mentioned 
above, from the bottom of the figure towards the top. In general, however, most of the 
surface near the notch was flat with no significant features or texture. This flat surface 
leads to the conclusion that significant rubbing took place. The rubbing “polished” the 
surface and removed all features that would have indicated the mode of fracture, e.g. 
fatigue striations, dimples, etc. The mbbed surface was visible over approximately 
80% of the surface when moving away from the notch toward the ridge of the fracture. 
Point B in Figure 6.3 is the approximate location where Figure 6.6 was taken. Figure 
6.6 is an example of the typical surface appearance in the rubbed region.^ The 


* Figure B. 1 in Appendix B also contains magnified views of the striated region. 

Figure B.2 in Appendix B is another picture of the typical surface appearance in the rubbed region. 
Point B in Figure 6.4 is the approximate location where Figure B.2 was taken. 
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polishing might have resulted from rubbing of the crack faces while the gear was in 
operation or rubbing against a part of the gearbox after fracturing away from the 
pinion. The extent and uniformity of the flat, polished surfaces support the former 
hypothesis. 



Figure 6.5: Fatigue striations near EDM notch on tooth #1 1 at 30x (a) and at 307x (b). 

A transition from the flat, polished area to one with some texture combined 
with flattened areas was observed further from the notch near the ridge (point C in 
Figures 6.3 and 6.4). In Figures 6.3 and 6.4, this combination, or partially rubbed, 
type of surface was found along the transition line from the darker region (flat, 
polished area) to the lighter region of the upper left comer. Recall the light region in 
both figures is near the toe end of the tooth. Figure 6.7, taken from region C in Figure 
6.3, shows clearly the features of the partially mbbed surface. The appearance of the 
raised areas is as if they have been flattened, while the lower lying regions have 
morphology indicative of fatigue. However, no well developed fatigue striations are 
observed. 
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Figure 6.6: Typical picture of flat, polished area on tooth #5 (410x). Photograph was 

taken near location B in Figure 6.3. 



Figure 6.7: Typical picture of partially rubbed surface (695x). Photograph was taken 

from location C in Figure 6.3. 

The lighter region in the upper left comers of Figures 6.3 and 6.4 (point D) 
shows little to no signs of rubbing. The surface also shows no obvious signs of 
fracture mode, e.g. intergranular fracture, ductile mpture, dimpling. Although no 
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striations were present on the surface, the fact that no other obvious signs of failure 
mode were observed, lead to the conclusion that the crack propagated by fatigue in 
this region under an applied load range which was inadequate to produce striations. 
The lack of rubbing also suggests that the fracture surfaces were created in the later 
stages of crack growth. 

The combination type of surface was also found over approximately 90% of 
the surface on the load free side of the tooth. Figure 6.8, of a partially rubbed surface, 
was taken from location A in Figure 6.9. Region B in Figure 6.9 is the tooth surface 
on the load free side. Therefore, point A is approximately 0.75 mm from where the 
crack ended on the tooth surface. In addition, there are fatigue striations evident in 
Figure 6.8. Because this figure is from the convex side of the tooth, the crack growth 
direction was from the top of the figure to the bottom. This combination of evidence 
leads to the conclusion that the crack continued to grow in fatigue mode along the 
convex side of the tooth. 

A light band can be seen in Region C of Figure 6.9. The darkened region 
separating region B and C is assumed to be oxidation of the fracture surface. Recall 
that the fatigue striations in Figure 6.8 are from location A. The surface in Region C 
shows obvious signs of ductile rupture. Figure 6.10. This observation is encouraging 
because it demonstrates that the material is capable of failing by ductile fracture, and 
the areas where this type of fracture occurred should be obvious and visible under the 
SEM. This result also leads to the conclusion that the primary mode of crack growth 
on the concave and convex sides of the tooth was fatigue. 



Figure 6.8; Picture of partially rubbed surface with fatigue striations on load free side 
(825x). Photograph was taken at location A in Figure 6.9. 
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Figure 6.9: Low magnification (31.4x) view of tooth #5’s fracture surface on convex 

side. 


A third tooth (#9) was also observed with the SEM. All of the features 
observed on teeth #5 and #11, with the exception of the ductile fracture area, were 
observed on tooth #9. No additional features could be seen. It is concluded that the 
observations made of teeth #5 and #11 are good representations of the crack patterns 
on all of the fractured teeth. 



Figure 6.10: Magnified view of ductile rupture at location C in Figure 6.9 (1670x). 
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Figure 6.11 summarizes the surface appearance on the loaded and load free 
sides of the fractured pinion teeth. A scenario of crack growth progression is 
developed based on these observations. The fatigue crack growth initiates from the 
EDM notch. The growth continues into the rim and at a larger rate towards the toe 
than the heel since it is assumed that the rubbed areas are the older surfaces. Once the 
crack reaches the ridge, the crack continues to grow toward the toe end. Figure 6.12 is 
a sketch of this scenario on the loaded side. The numbers in the sketch correspond to 
the progression of the crack front. When the cra^ reaWes the tooth surface at the toe 
end, the extent of crack growth has dramatically changed the stress distribution in the 
remaining ligament. Consequently, the crack front turns toward the fillet on the 
convex side, and progresses by fatigue along the convex side. When the crack front 
becomes sufficiently close to the root of the convex side, ductile rupture occurs in the 
remaining ligament. After this, any additignal load on the tooth causes the torsional 
tearing of the ligament on the heel end. Figure 6.2 sketches the crack growth through 
the tooth width. This sketch is applicable to cross-sections from the toe end to 
approximately the middle of the tooth length. 



Figure 6.1 1 : Sketch of loaded and load free sides of a pinion tooth’s fracture surface 
appearance along the length. Orientation is consistent with SEM pictures. 
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Figure 6.12: Sketch of crack propagation scenario on loaded side devised from 

fracture surfaces. 


6.4 Chapter Summary 

This chapter was devoted to presenting data from an OH-58 spiral bevel pinion 
test. The test was conducted by NASA/GRC. EDM notches were introduced into the 
root of nine of the pinion teeth to serve as starter cracks for fatigue crack growth. 
Limited observations of the crack growth during the test were made, and, as a result, 
the fracture surfaces were observed with a SEM. 

Overall, the microscopy identified fatigue crack growth regions and regions of 
ductile rupture successfully. In addition, the crack face morphology showed 
significant signs of rubbing, which had “polished” the surface. This polishing 
removed any discemable fracture surface features on the majority of the surfaces. The 
signs of fatigue on the loaded and load free sides of the fracture surface indicated the 
majority of the crack growth was attributed to fatigue. At the ridge near the toe end, 
the surface showed little to no signs of rubbing. This observation suggested that the 
surface was created in the latter stages of crack growth. It was inferred from the 
jagged and tom appearance of the fracture surface near the heel that this region was 
the last remaining ligament of the tooth after rupture occurred in the root of the convex 
side. 

Due to the dearth of well-developed fatigue striations on the fracture surfaces, 
no observations were made on the crack growth rates. In addition, the large amounts 
of rubbing removed all indications of crack front shape during propagation. 
Nevertheless, a scenario of crack propagation was devised. The next chapter 
compares these test results to the simulations results from Chapter 5. 
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CHAPTER SEVEN; 

DISCUSSION AND SENSITIVITY STUDIES 


7.1 Introduction 

In Chapter 5, fatigue crack growth in the OH-58 spiral bevel pinion under 
moving, non-proportional loads was predicted. Chapter 6 presented experimental 
fatigue crack data from a tested OH-58 spiral bevel pinion. The present chapter 
compares these two sets of results to evaluate the success of the predictions and 
investigates the sensitivity of the prediction results to variations in the methods 
assumptions. 

In Section 7.2, the crack growth simulations are compared to the experimental 
results of a tested pinion; the fatigue lives and crack trajectories are evaluated to 
determine the accuracy of the prediction method. Sensitivity studies are conducted in 
Section 7.3 to explore variations in tooth contact position and magnitude and the 
sensitivity of the crack-closure-based fatigue crack growth rate models to variations in 
the model parameters. 

Crack growth predictions from the moving load analyses (Section 5.5) are 
compared to crack growth results from analyses that consider only highest point of 
single tooth contact (HPSTC) loading in Section 7.4. HPSTC loading has been 
commonly adopted in past research because it is a more simplified approach than the 
moving load. When using HPSTC loading, existing fatigue crack growth theories can 
be implemented since there is a single load location and proportional loading. The 
two loading methods’ results are compared to evaluate the need for the moving, non- 
proportional load method; the least computationally intensive model and method 
which produces reasonable crack growth results is the most practical for a gear 
designer. 

7.2 Comparisons of Crack Growth Results 

Figure 7.1 shows the predicted and experimental trajectories on the tooth 
surface and through the cross section of a pinion tooth. The predicted results are the 
same as those reported in Section 5.5. The experimental trajectories are approximated 
from measurements and photographs of the failure surfaces of tooth #11. As 
discussed in Chapter 6, the failure associated with tooth #1 1 is representative of all the 
failures in the tested pinion. In addition, the size of the initial flaw in the predictions 
was taken from the dimensions of the EDM notch in tooth #11. 

A ridge is not observed in the predicted crack path through the cross section 
because the simulations were stopped prematurely. The simulations were halted 
because the trajectory along the tooth surface on the toe end varied significantly from 
the experiment. As a result, it was concluded that the simulations were not completely 
accurate and further propagation of the crack would lead to no additional insights. An 
explanation for this discrepancy might be that the simulation did not properly account 
for changes in the load shape and location as the crack grew. This changing load 
scenario will be investigated further in Section 7.3.3. 
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Experimental 

\ V. Predicted / 


b) Cross section at location of midpoint of initial crack 
Figure 7.1 : Predicted and experimental crack trajectories. 


Qualitative comparisons of the analyses’ results to the tested pinion results can 
be made. The numerical analyses predict tooth failure, which is concordant with the 
experiment. In addition, in both the fractography study and numerical analyses, the 
crack propagates more rapidly toward the toe than the heel. In both the test and 
simulation, a portion of the tooth at the heel end remains intact. Additionally, the 
simulations predict the crack propagating along a steeper trajectory into the gear rim in 
the middle of the tooth length than on the toe end of the length. This behavior is also 
observed in the tested pinion. 

The final predicted trajectory through the thickness of the tooth agrees very 
well with the initial path in the experiment. It is assumed that this path could lead to 
the formation of a ridge if the simulations were continued. The entire predicted crack 
trajectory, however, does not completely match the tested pinion. The simulations 
predict the toe end of the crack turning up the tooth height at a steep angle. This 
behavior is not seen in the tested pinion. One reason for the discrepancy could be that 
the loading conditions for the simulations were not identical to the test. The gear was 
tested at increasing torque levels over the 4.9 million cycles. The simulations, 
however, were performed under a constant torque level. The increase in torque should 
affect the tooth contact, which, in turn, will influence the crack trajectory. The 
influences of the torque level and contact location on crack trajectories are explored in 
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Section 7.3.3. The goal is to predict what operating conditions would be necessary to 
obtain the fracture path observed in the tested pinion. 

For the constant torque level of 3120 in-lb, approximately five hundred 
thousand cycles were predicted to propagate the crack through the tooth thickness. 
This value is the same order of magnitude as that which occurred in the tested pinion. 
The number of cycles to produce the tooth fractures in the test is smaller than 4.9 
million. The sensitivity of the fatigue life prediction to the values chosen for the 
model parameters is studied in Section 7.3.1. 

Alban’s condition number four for “classic tooth-bending fatigue” scenario 
(Section 2.5) is not captured in the numerical work. The magnitude of the applied 
loads during the cycle was kept constant during the crack propagation analyses. Figure 
7.2a. This type of loading scenario is considered load control and results in the SIFs 
increasing continuously as the crack grows. Figure 7.2b. Fatigue crack growth will 
occur at an increasing rate until the SIFs satisfy a fracture criterion, such as K j — K . 
Crack growth simulations under load control will predict a shorter number of cycles to 
failure than observed. This is because, in reality, when a cracked and uncracked tooth 
mesh, the cracked tooth will deflect a limited amount before it’s adjacent tooth picks 
up a portion of the load [Alban 1985]. The displacement of the crack faces reaches a 
maximum and will be roughly equal for every remaining load cycle. As a result, the 
rate of increase in the SIFs will decrease, and reach roughly a constant maximum 
every cycle. An idealization of this is shown schematically in Figure 7.2c. 
Propagating the crack under these conditions is considered displacement control. 
When the maximum SIF ceases to increase, the fatigue crack growth rate is relatively 
constant, and the number of cycles to failure increases. 



control, and c) Ki versus time, displacement control. 
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To demonstrate that the simulations are capable of predicting the turning that is 
necessary to predict the ridge in the fractured tooth, a large crack shape is assumed 
and inserted into the full pinion BEM model. The crack front coordinates are 
determined by the location of the ridge in a fractured tooth from the experiment. 
Figure 7.3 is a photograph of the fracture surface with the approximate location of the 
assumed crack front designated by the dashed line. Based on the SEM observations, 
the assumed crack has propagated along the root from the initial notch to the toe 
surface. Figure 7.4 is a picture of the BEM geometry model that illustrates the 
assumed initial crack trajectory on the tooth surface. Since the correct contact areas 
for a tooth that is flawed to this large of an extent are unknown, the HPSTC load step 
from the moving load analyses (load step eleven) is used. 



Figure 7.3: Assumed location of crack front (ridge). 



Figure 7.4: Tooth surface showing assumed shape of large crack (dashed line). 

A few cycles of crack growth are carried out using the method described in 
Section 3.2.3. This method assumes mode I dominant fatigue crack growth with 
static, proportional loading. The cycles of crack growth are necessary to demonstrate 
the direction the crack front progresses from its assumed location at the beginning of 
the ridge formation. Figure 7.5 shows the trajectory through the thickness of the tooth 
at approximately the middle of the tooth length. The initial trajectory into the rim is 
assumed to be flat, and the curving at the end of the crack length shows the formation 
of the ridge. This demonstration of the ridge formation based on an assumed crack 
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front from the SEM observation supports the crack growth scenario developed in 
Section 6.3.2. 


/ \ 



Figure 7.5; Crack trajectory through tooth thickness for assumed large crack. 

One could argue that the tested spiral bevel pinion failure did not exactly meet 
the classical failure conditions described in Section 2.5 since the entire length of the 
tooth did not fracture from the gear; a portion of the tooth at the heel remained intact. 
This suggests that the loading might have been biased toward the toe end of the tooth. 
Numerical analyses with shifted load locations are presented in Section 7.3.3 that give 
insight into the sensitivity of the crack trajectories to loading location. 

The predictions in Section 5.5 did not consider changes in the original contact 
locations during propagation. The increasing tooth deflections as the crack grows 
might cause the original contact locations to shift and the distribution of load and the 
size of the contact ellipses to change. A three dimensional, contact mechanics, and 
fracture mechanics simulation of the rolling process between two mating gears is 
necessary to capture the load redistribution effects fully. This type of analysis is not in 
the scope of this work. 

It is impossible to determine the exact amount of rubbing between the crack 
faces based on the BEM analyses. In Section 6.3.2 it was concluded that the surfaces 
with greater amounts of rubbing were formed in the earlier stages of crack growth; 
since these surfaces were older, the features of the surfaces had more time to rub 
away. The kinematics of the geometry and loading is another explanation for the 
varying amounts of rubbing observed on the fracture surfaces. Rather than attribute 
the varying degrees of rubbing to time, it could be attributed to the magnitude of 
contact between the fracture surfaces. The loading might deflect the tooth in a manner 
that does not allow the ridge’s fracture surfaces to rub, but does create large contact 
forces between the crack faces near the notch. A three dimensional analysis modeling 
contact between the crack faces with accurate loading conditions on the tooth surface 
is necessary to determine the true cause of the rubbing. 

7.3 Sensitivity Studies 

These studies are performed to gain insight into the sensitivity of predicted 
crack growth rates and predicted crack trajectories to growth rate model assumptions, 
load magnitude, and load location. They are also conducted to investigate possible 
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causes for the discrepancies between the predictions and the experimental results. The 
fatigue crack growth rate model parameters. Section 7.3.1, crack closure model 
parameters, Section 7.3.2, and the contact position and magnitude. Section 7.3.3, are 
researched further. 

7.3.1 Fatigue Crack Growth Rate Model Parameters 

Limited fatigue crack growth rate data is available in the literature for AISI 
9310 steel. The predictions in Section 5.5 used values for the crack growth rate model 
parameters, n and C, taken from a curve fit to the intrinsic fatigue crack growth rate 
data (no closure) for AISI 9310 tested in 250° oil. These values were 3.36 and 
6.19e-20 (in/cycle)/(psi*in° ^)", respectively. The range of values for n and C from the 
literature is reported in Table 7.1. Au et al.'s data are not from intrinsic fatigue crack 
growth rate curves. Their data are from fatigue crack growth tests with ^ = 0.05 . 
The other three sets of model parameters have been normalized to an intrinsic fatigue 
crack growth rate curve. 


Table 7.1: F atigue crack growth rate model parameters. 


Source 

n 

C 

[(in/cycle)/(psi*in°'^)]'' 

[in/cycle] 

Cycles 

/inch* 

Fatigue 

Life* 

Forman et al. 
[1984] 

1.63 

1.08e-13 

9.30e-6 

107,527 

675,838 

Au et al, 
[1981]" 

2.56 

2.72e-17 

2.03e-6 

492,611 

3,096,201 

Air test 

[Proprietary 1998] 

3.63 

5.49e-21 

1.54e-5 

64,935 

408,145 

Oil test 

[Proprietary 1998] 

3.36 

6.19e-20 

1.23e-5 

81,301 

511,000 


The fatigue life estimates using each set of parameters in Table 7.1 assume that 
the number of cycles that each source would predict for the gear’s fatigue life is 
roughly proportional to the ratio of the oil test’s cycles/inch to each set’s cycles/inch. 
Forman et al., Au et al., and the air test data each predict a fatigue life of 675,838 
cycles (32% increase), 3,096,201 cycles (506% increase), and 408,135 cycles (20% 
decrease), respectively. 

Au et al.'s combination of n and C predicts the smallest growth rate and 
therefore the longest fatigue life. The benefits and conservatism of considering crack 
closure in the predictions is demonstrated by comparing the predictions using the 
intrinsic parameters to the predictions with Au et al.'s parameters. Figure 7.6 contains 
the fatigue crack growth rate curves for the various sets of n and C. The curves are 
generated from the data in Table 7.1. 


Calculations are based on an assumed value for AK - 18,000 psi*in” \ 

" This data was taken from a fit to fatigue crack growth rate data for non-carburized AISI 9310 tested 
in wet air, at a loading frequency of 1 .0 Hz, and R = 0.05 . The parameters are not from intrinsic 
fatigue crack growth rate data. It is should be noted that when parameters from the air test at R = 0.05 
are used the calculated growth rate is 2.02e-6 in/cycle. 
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Figure 7.6; Fatigue crack growth rate curves for the sets of model parameters in 

Table 7.1. 

An exact evaluation of which set of material constants is most accurate is 
nearly impossible since many parameters of the pinion test are unknown. The 
calculations presented in this section demonstrate a trend that the fatigue life 
calculations will be more accurate when material constants from intrinsic fatigue crack 
growth rate curves are used. 

7.3.2 Crack Closure Model Parameters 

For the simulation results presented in Section 5.5, values for the crack closure 
model parameters, jcand /?, were assumed in order to calculate the fatigue crack 
growth rates. This section investigates the validity and sensitivity of the results to the 
assumed values. 


Crack Growth Rate Sensitivity to K 

K incorporates three dimensional effects into the crack growth rate 
calculations. Newman specifies that jc varies between one and three for plane stress to 
plane strain, respectively. For the predictions reported in Section 5.5, K was equal to 
three. However, for extremely shallow cracks or portions of the crack front near the 
free surface, a value of k equal to one might represent the crack conditions more 
accurately. 

One method to evaluate the crack tip conditions is to compare the size of the 
crack tip plastic zone to the crack’s geometry. The plastic zone is larger in plane 
stress than in plane strain. An approximation of the plastic zone size, fp, is 




J_ 

3;r 


(T. 


V 


(7.1) 
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Using load step eleven’s SIF results from the initial crack (Figure 5.8), which 
are the largest mode I SIFs during the load cycle, the plastic zone size along the crack 
front ranges from 1.44x10 '* inches to 1.07x10'^ inches. These dimensions are only 
0.29% and 2.14%, respectively, of the initial crack depth. It is concluded, therefore, 
that the plane strain assumption along the entire crack front most accurately represents 
the conditions in the real gear. 

Crack Growth Rate Sensitivity to ^ 

S„,a.x, the far field applied stress, is a function of P,c (half of the crack length), 
and K/. ^ is a dimensionless quantity that considers geometry effects when relating Ki 
to the applied stress. Values of /? from handbook solutions can vary from one half to 
two [Murakami 1987]. Since the gear geometry is complex and unlike any handbook 
solution, a value of /3 = 1 was selected. An alternate approach could have been to use 

a known P factor for a similar, simplified geometry. This alternate approach will now 
be investigated and growth rates between the two methods will be compared. 



Figure 7.7: Finite thickness plate with a semi-elliptical surface crack subjected to 

mode I uniform stress. 

The initial crack in the gear is approximated by a finite thickness plate with a 
semi-elliptical surface flaw subjected to mode I uniform tensile stress, <7, Figure 7.7. 
The magnitude of Ki varies along the crack front. and occur at the surface 

and midpoint of the crack front, respectively. They are given by Broek [1986] as: 

= and (7.2) 
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, ^ 3;r K 

where o 1 — 

8 8 




) 


From Equation (7.2), the maximum and minimum expressions for are 

„ 1.12 1.12 

/3„„= — and/3.„,=— 


(7.3) 


Based on the initial crack geometry, P„,a,y and /3„„„ are 0.783 and 0.701, respectively. 
Kop and da are recalculated using these values for (3, the SIFs from the initial notch 
analyses (Figure 5.8), and the same model parameters as were used in Section 5.5. 



Crack front position 
(Orientation; heel to toe) 

Figure 7.8: Change in K^p and da as functions of P„a.\ and Pmin with respect to original 

calculations with P = 0. 

Figure 7.8 shows the percent change of Kop and da with respect to the original 
calculations. The data show that the largest percent difference is 2.2%. As the crack 
grows the ratio of o to c will become smaller since the tooth length is longer than the 
tooth width. This will increase (p and, therefore, increase P„,ax and decrease 
However, as the crack grows, this closed form solution for K/ in the gear is no longer 
valid since the crack’s geometry changes dramatically. Therefore, no further 
conclusions can be stated on the effect of P on crack growth calculations. 

7.3.3 Loading Assumptions 

The intent of this study is to determine how the crack trajectory changes under 
different contact conditions. One motivation for this is that the simulation and 
experiment’s crack trajectories on the toe end do not match. The tested pinion’s crack 
mouth remained relatively flat along the root until it reached the end of the tooth 
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length at the toe (Figure 6.1b). The trajectory in the simulation turned, out of the root, 
up the tooth height and eventually reached the top land (Figure 5.14a). It is 
hypothesized that the differences may be attributed to inconsistencies between the 
contact conditions (loading conditions). 

The inconsistencies could result from misalignment during the test or 
inaccurate representation in the simulation of the actual contact areas in the test. 
Glodez et al.'s [1998] experimental work with spur gears supports this hypothesis. 
They considered two load cases: i) loading along the entire length of the tooth, and ii) 
loading along one half of the length. With load case ii) the crack in the unloaded 
portion of the tooth length turned out of the root and grew up the tooth height. On the 
other hand, the crack in load case i) remained flat along the entire length of the tooth 
root. The goal of the remainder of this section is to investigate whether shifted loads 
have the same influence on crack trajectories in spiral bevel gears as Glodez et al. 
showed in spur gears. 

Load Location 

Two shifted load scenarios are investigated. For both scenarios, the cracked 
BEM model from propagation step number five is analyzed under the shifted contact. 
This model was chosen because the crack trajectory began turning sharply from this 
step onward in the preliminary analyses. The contact areas are shifted approximately 
+0.3 inches along the tooth length. 

The crack trajectory for the shifted contact areas is calculated using the non- 
proportional load method described in Section 5.4.2. The load cycle is approximated 
by the discrete load steps one, five, and eleven. Table 7.2 sketches the predicted 
trajectories on the tooth surface for the original and two shifted analyses. The mode I 
and n SIFs from the shifted loading scenarios are given in Appendix C. 

When the contact is central, the crack turns up the tooth height on both ends. 
The trajectory “wraps around” the contact location. However, when the contact is 
shifted toward the heel (toeX the tendency for the crack to kink up on the heel end (toe 
end) is suppressed. This is most clearly seen when comparing the central and toe 
contact location trajectories. This result is consistent with Glodez et aids 
observations. As a result, it is assumed that, if the fatigue crack growth simulations 
were carried out further with the shifted contact locations, a flatter trajectory that 
maintains a path very near the root under the contact location will result. 

The discrepancy of the toe end trajectory between the test and simulation is 
explained by the fact that, in the test, the contact was closer to the toe end. The shifted 
contact could have resulted from increasing deflections of the tooth. As the crack 
grew, the tooth’s stiffness decreased, and the load could have been redistributed along 
the tooth length. The subtleties of the redistribution and its effect on crack trajectories 
can only be modeled accurately with a three dimensional contact analysis between the 
mating gears in conjunction with a fracture mechanics simulation. 
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Load Magnitude 

The tested spiral bevel pinion was run at varying levels of input torque detailed 
in Table 6.1. However, the simulation results reported in Section 5.5 assumed contact 
areas and load magnitudes produced by 3120 in-lb torque (100% design load). The 
goal of the current study is to identify the influences of the increased torque levels on 
crack trajectories. 

The SIF distributions and trajectories under larger torque levels of 3874 in-lb 
and 4649 in-lb (125% and 150% design load, respectively) are explored. From 
Hertzian contact theory, it is known that the lengths of the contact ellipses’ axes are 
proportional to the cube root of the applied load (Equation (2.1)). Consequently, the 
lengths of the major and minor axes incr ea se by 7.72% (125% design load torque) and 
14.47% (150% design load torque) under the larger loads. It is assumed that the mean 
contact points (center of the ellipses) are the same as the points defined for the 100% 
design load. Similar to the shifted load analyses, the crack from propagation step five 
in the moving load simulations (Section 5.5) is selected to analyze under the larger 
torque levels. Figure 7.9 shows the locations of the crack and of the contact ellipses 
defined for 125% design load. 
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Figure 7.9; Geometry model with crack showing contact areas one and eleven defined 

for 125% design torque. 


The mode I and II SIFs from load steps one and eleven are presented in Figures 
7.10 and 7.1 1, respectively. The mode I SIFs do not increase linearly with the larger 
loads; the 125% load has a larger effect than the 150% load. The smaller spread 
between the curves produced by load step one at the toe end is most likely explained 
by the fact that the load is not over this portion of the crack. In contrast, the SIFs 
increase uniformly along the entire crack front for load step eleven; the major axis of 
contact ellipse eleven is larger than the length of the crack mouth, and the ellipse is 
located directly above the crack. On the other hand, for load step one the influence of 
the increased load on the mode II SIFs is opposite. There is a larger spread in the 
curves over the portion of the crack with no load above it (toe end). 

The ratio KulKi is important because it determines the crack trajectory angle 
and the amount of rubbing between the crack faces. The larger the ratio is, the larger 
the kink angle will be and the greater the amount of rubbing. Figure 7.12 contains 
these ratios produced by the two load locations and all three load magnitudes. The 
curves demonstrate that the ratio of Ku to Ki increases as the magnitude of load and 
size of the contact area increases. This implies that Ku is more sensitive to the 
changes in the torque level than Ki. This result supports the fractography 
observations. A large percentage of the fracture surface displayed signs of significant 
amounts of rubbing between the crack faces. The SEF ratios from the initial crack 
propagation analyses were not necessarily large enough to support the extent of 
rubbing observed. However, it appears that the increased torque levels will increase 
the amount of rubbing between the crack faces. 

The kink angles calculated by the maximum principal stress theory for the 
various load locations and magnitudes are given in Figure 7.13. The largest absolute 
change in angle is 9.7° and 6.4° for load step one and eleven, respectively. 
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Figure 7.10; Kj distribution for load step one (a) and load step eleven (b). 
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Figure 7.12: KulKi distribution for load step one (a) and load step eleven (b). 


NAS A/CR— 2000-2 1 0062 


95 





Figure 7.13: Kink angle distribution based on the maximum principal stress theory for 
load step one (a) and load step eleven (b). 


7.4 Highest Point of Single Tooth Contact (HPSTC) Analysis 

Comparison studies to determine the smallest model that accurately represents 
the gear’s operating conditions were performed when developing the BEM model. 
These results were reported in Section 5.2. Similar comparisons are now made 
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between the moving load method and a simplified loading method. Again, the 
assumption is that the moving load method is most accurate. The simplified method 
assumes a cyclic load at the HPSTC on the pinion tooth. The HPSTC is taken as 
contact ellipse number eleven from the discretized moving load; load step eleven is the 
final step of single tooth contact in the discretized load data. The magnitude of the 
load is defined as 100% design load. The model parameters and material properties 
from the moving load analyses are used in the HPSTC predictions. 

The initial crack location and geometry are the same as those from the moving 
load analyses. The method to propagate the crack under the HPSTC is described in 
Section 3.2.3. The method assumes proportional loading. It is assumed that the 
HPSTC produces Kj,na, and that R is zero. The direction of growth is determined by 
the maximum principal stress theory using the ratio of Kn to Kt from the HPSTC 
loading. The extensions for the discrete crack front points are calculated with Paris’ 
model modified to account for crack closure. Figure 7.14 is a comparison of the crack 
trajectories from the moving and HPSTC load methods. Roughly 190,000 cycles have 
occurred. The cross section view is taken at the approximate location along the tooth 
length of the initial crack front’s midpoint. 



a) Tooth surface 



b) Cross section of tooth at midpoint of 
initial crack 

Figure 7.14: Comparison of crack trajectories from moving load and HPSTC load 
(fixed location) methods after ~ 190,000 cycles. 
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The midpoint of the crack front is deeper in the HPSTC analyses after 190,000 
thousand cycles. From Figure 7.14b, it appears that the moving load analysis 
trajectory will produce rim failure. Figure 5.14, however, shows that the crack turns 
when the predictions are continued. Therefore, both the static and moving load 
method predict tooth failure. The slope of the trajectory into the rim in the moving 
load prediction matches more closely the observed trajectories in the tested pinion. 
This comparison is purely qualitative. 

Several obvious differences between the trajectories predicted by the two 
methods can be observed. As seen in Figure 7.14a, the HPSTC method predicts a 
larger kink at the , heel end; the moving load method predicts a larger kink at the toe 
end. Considering the location of the HPSTC load, this result is consistent with the 
shifted load analyses of Section 7.3.3. 

One may conclude from Figure 7.14b that the HPSTC method predicts a larger 
crack face area since the cross section view of the crack is deeper, yet the lengths of 
the cracks on the tooth surface are roughly equal. Figure 7.15, in general, supports 
this conclusion. 


0.25 : 



N [cycles] 

Figure 7.15: Crack area versus number of load cycles for HPSTC and moving load 

prediction methods. 

In summary, the HPSTC analyses predict the same failure mode as the moving 
load analyses. The crack trajectory and fatigue life calculations vary between the two 
methods. Since no experimental fatigue life data exists, the accuracy of one methods 
fatigue life prediction over the other methods can not be evaluated. The moving load 
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predictions match the experimental trajectory into the rim and through the cross 
section of the tooth better than the HPSTC prediction. Since the trajectory into or 
through the rim is what determines tooth failure or rim failure, it is concluded that the 
moving load method is necessary to capture that result most accurately. All of the 
trajectories on the tooth surface at the heel end, however, are in reasonable agreement. 
Nonetheless, a distinct advantage of the HPSTC method is the significant decrease in 
computational time to perform the crack propagation predictions since only one load 
case needs to be analyzed. 

7.5 Chapter Summary 

The results from a fatigue crack growth simulation in a spiral bevel pinion 
were compared to crack growth observations in a tested pinion. The comparisons are 
summarized as follows: 

• The simulations predicted a reasonable fatigue life with respect to the test data. 

• The original trajectory predictions failed to capture detailed aspects of the 
observed fracture surfaces in the test. It was determined that the simulated loading 
on the tooth probably modeled the tooth contact in the test incorrectly. The tooth 
contact information used in the predictions assumed perfect alignment between the 
pinion and the gear and that the gears were not flawed. Some explanations for the 
differences in contact between the test and theory were determined to be: 

1. Change in contact location in the test as the crack grew and the tooth 
became more flexible. 

2. Differences in the magnitude of loading. 

3. Crack growth under load control (simulation) versus displacement control 
(test). 

4. Misalignment between the gear and pinion in the test. 

• Additional simulations demonstrated the capability to predict the crack trajectory 
observed in the test. A large initial crack, which was assumed to approximate the 
location of the crack front just prior to the formation of the ridge, was used and the 
crack was propagated through a series of steps. 

Sensitivity studies were conducted to determine how changes in some of the 
crack growth method’s assumptions would modify the predictions. The studies 
determined that: 

• The fatigue crack growth rate model parameters used in the initial prediction 
yielded conservative results. 

• The crack front condition is best described as plane strain. 

• A reasonable approximation of the dimensionless quantity j3, which incorporates 
geometry effects when calculating SEFs, is a value of 1 . 

• The trajectory observed in the tested pinion would result from a contact biased 
toward the toe end. 

• The increased torque levels might explain the significant amounts of rubbing seen 
on the fracture surfaces of the tested pinion. 
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A simplified loading method that assumes a cyclic load at the HPSTC on the 
pinion tooth during meshing was investigated. The failure mode predicted by this 
method was the same as the moving load predictions. However, the crack trajectory 
and fatigue life calculations varied between the two methods. The HPSTC method is 
advantageous because it significantly reduces the computational time. However, upon 
comparison of the results from the two methods to experimental results, it is 
concluded that the moving load method’s trajectories are more accurate. 

In summary, insights into the intricacies of modeling fatigue crack growth in 
three dimensions were gained. Preliminary steps toward accurately modeling crack 
growth in complicated three dimensional objects such as spiral bevel gears were 
completed successfully. To improve the accuracy of the simulations, the change in 
contact between spiral bevel gear teeth during operation as a crack evolves is needed. 
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CHAPTER EIGHT: 
CONCLUDING REMARKS 


8.1 Accomplishments and Significance of Thesis 

This thesis investigated computationally modeling fatigue crack growth in 
spiral bevel gears. Predicting crack growth is significant in the context of gear design 
because a crack’s trajectory determines whether the failure will be benign or 
catastrophic. Having the capability to predict crack growth in gears allows a designer 
to prevent catastrophic failures. Prior to this thesis, numerical methods had been 
limited to modeling cracking in gears with simpler geometry, such as spur gears. Spur 
gear geometry permits the use of two dimensional analyses. However, spiral bevel 
gears require a three dimensional model of the geometry and cracks. Three 
dimensional models are much more complicated to create, require greater computing 
power because of the significant increase in degrees of freedom, and no closed form 
solutions exist to predict the growth of arbitrary three dimensional cracks. Prior to this 
thesis, few predictions of crack growth in spiral bevel gears had been performed. 
Accurately modeling three dimensional fatigue crack trajectories in a spiral bevel 
pinion required expanding the state-of-the-art capabilities and theories for predicting 
fatigue crack growth rates and crack trajectories. 

The geometry of a spiral bevel pinion from the transmission system of the U.S. 
Army’s OH-58 Kiowa Helicopter was used for demonstrative purposes. A BEM 
model of the pinion was developed using a computer program developed by 
NASA/GRC that calculates the surface coordinates of a spiral bevel gear tooth. Their 
tooth contact analysis program was also used to determine the location, orientation, 
and magnitude of contact between the pinion and its mating gear. The contact was 
represented by discrete traction patches on the gear tooth. 

LEFM theories were combined with the BEM to accomplish the crack growth 
predictions. The simulations were based on computational fracture mechanics 
software developed by the Cornell Fracture Group, which allow for arbitrarily shaped, 
three dimensional curved crack fronts and crack trajectories. The crack trajectories 
were determined by a Paris model, modified to incorporate crack closure, to calculate 
fatigue crack growth rates in conjunction with the maximum principal stress theory to 
calculate kink angles. 

In operation, the load on a gear tooth varies over time in location and 
magnitude. This moving load effect was incorporated into the propagation method. 
Only loads normal to a gear tooth’s surface were considered. It was discovered that 
the moving normal load produces a non-proportional load history in the tooth root. 
Proposed prediction methods for fatigue crack growth under non -proportional loads in 
the literature were determined to be insufficient for the spiral bevel gear model. As a 
result, a method to predict three dimensional fatigue crack growth under non- 
proportional loading was developed. The method incrementally advanced the crack 
front for a series of discrete load steps throughout one load cycle. A number of load 
cycles were then specified, and the crack was advanced an amount based on the 
number of specified load cycles and the calculated trajectory from the single load 
cycle; the process was then repeated. Some aspects of the final crack trajectory 
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predicted by this moving load method differed from a failure in a tested pinion; 
however, the method succeeded in predicting a fatigue life comparable to the 
experimental data. 

Other issues related to modeling crack growth in a gear were also investigated. 
For example, the effect of shifting the load location along a tooth’s length on the crack 
trajectories was confirmed. For a crack that has initiated in the tooth’s root, when the 
load location is directly above the crack, the crack trajectory will remain very close to 
the root. Additionally, the effect of compressive loads on fatigue crack growth rates in 
AISI 9310 steel was examined. This examination is significant because a principal 
focus of current gear design is to minimize a gear’s weight. Reducing the amount of 
material in the gear may increase the magnitude of the compressive stresses in a gear 
tooth’s root, which could influence crack growth rates. It was discovered that the 
compressive portion of a load cycle did not significantly modify the rates when crack 
closure was incorporated into Paris’ model to calculate fatigue crack growth rates. As 
a result, the BEM/LEFM analyses of a spiral bevel pinion were carried out ignoring 
the compressive portions of the loading history. 

The predictions from the moving load crack propagation method were 
compared to predictions when only HPSTC was considered. HPSTC is a more 
simplified approach and has been commonly used in past research when numerically 
analyzing crack propagation in gears. The HPSTC method utilized existing fatigue 
crack growth theories since there was a single load location and proportional loading. 
The analyses in this thesis with the two loading methods predicted different fatigue 
lives and crack trajectories. The lack of experimental fatigue crack growth rate data 
hindered an evaluation of the crack growth rates predicted by the two methods. The 
moving load method’s crack trajectory predictions agreed more closely to the tested 
pinion failures. Crack trajectories are of primary importance to predict the failure 
mode. 

The dearth of fatigue crack growth rate data and crack front shape information 
from tooth failures in a tested spiral bevel pinion motivated SEM observations of the 
fracture surfaces. A crack growth scenario was devised from the observations. In 
addition, the observations suggested that the failure mechanism along the majority of 
the surface was fatigue. This result supported the use of the numerical simulations to 
predict fatigue crack growth trajectories in the gear. 

As this thesis was a first attempt at predicting fatigue crack growth in spiral 
bevel gears, certain limitations were encountered. The limitations can be summarized 
as follows: 

• A scarcity of experimental data prohibited validations of calculated crack growth 
rates, fatigue life predictions, and crack front shape evolution. 

• The effect of tooth deflections on the contact area between mating gear teeth was 
not modeled. Capturing this effect will increase the accuracy of the model since 
crack trajectories are strongly determined by the load locations. 

• It is anticipated that the deflections of a cracked spiral bevel gear tooth will be 
limited by the adjacent tooth picking up the load. The magnitude of this 
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maximum deflection is unknown. As a result, the simulations did not model this 
behavior. 

• The method to predict crack trajectories under the moving load was time 
consuming since every discrete load step had to be analyzed for a full pinion 
model. For the method to be more practical to a gear designer, improvements in 
the computational time and implementation would have to be made. 

In summary, the achievements of this thesis brought the current gear design 
approach closer to developing a numerical tool to evaluate safety aspects of gear 
geometries. Critical areas that must be understood in greater detail prior to predicting 
more accurate crack trajectories and growth rates in three dimensions were identified. 
By incorporating non-proportional loading, complex three dimensional geometry, and 
arbitrary three dimensional crack models, the work extended the current capabilities 
and theories for predicting fatigue crack growth rates and crack trajectories. 

8.2 Recommendations for Future Research 

The accomplishments of this thesis produced many new questions and issues 
related to simulating fatigue crack growth in spiral bevel gears. Future research 
focused on the following areas will further assist the development of numerical design 
tools to evaluate a gear design’s safety. 

• A testing program centered on fatigue crack growth rates and trajectories from 
non-proportional loads should be carried out. A more fundamental understanding 
of crack behavior in this type of loading environment is needed to enhance the 
accuracy of fatigue crack growth rate models. To gain insights into the 
fundamental behavior, the tests do not have to be conducted with spiral bevel 
gears. To judge the correctness of the predictions in this thesis, however, more 
tests on spiral bevel gears that record crack growth rates and trajectories over a 
time period are necessary. Such test data are essential to confirm the accuracy of 
the simulations and to evaluate the proposed method for predicting crack growth 
under the moving tooth load. 

• A more detailed understanding of the contact between a cracked gear tooth mating 
with an uncracked gear tooth is required. The load redistribution effects in this 
scenario could be studied either experimentally or numerically. Numerical studies 
would require a fully three dimensional, LEFM, and contact analysis of the rolling 
process between two mating gears. Once the redistribution effects are captured, 
the findings could be applied to crack growth simulations. The expectation is that 
the simulations would model more accurately the observed behavior in real 
failures. 

• The analyses in this thesis considered only loads normal to the tooth surface. 
However, a gear tooth is subjected to a variety of loads in operation. For example, 
the contact between the mating gears also produces frictional forces along with the 
normal loads over the contact ellipse area. Gears in rotorcraft applications operate 
at elevated temperature, and, therefore, thermal effects might be included. In 
addition, the rotation of the gear produces centrifugal forces. Dynamic loads are 
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produced from the contact between the mating gears in conjunction with the high 
loading frequency. Several researchers have also modeled the residual stresses in 
a gear, which result from the difference in hardness between the case and core, in 
numerical analyses by thermal loads. Parametric studies to determine the 
significance of all these loading variables on predicted fatigue crack growth in a 
spiral bevel gear is imperative to the development of numerical gear design tools. 

• Further work to enhance the speed of the numerical analyses and increase the SIF 
accuracy would be beneficial. One method to improve the accuracy is to use the 
FEM. However, meshing three dimensional volumes with cracks introduces 
additional difficulties. Research is currently underway, and should be continued, 
to overcome these meshing difficulties. The computational time will continue to 
decrease as computer technology rapidly advances. The overall objective of 
developing a practical and accurate numerical design tool for any type of gear is 
foreseeable in the near future. 
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Figure A.l: Mode II SIF distribution for load steps one through eleven. 
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Figure A. 2: Mode III SIF distribution for load steps one through eleven. 
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APPENDIX B 



Figure B.l: Fatigue striations in Figure 6.5 at 360x (a) and l,000x (b). 

An additional point of interest is captured in Figure B.3. The figure is of the 
tooth surface. The horizontal lines are the grooves from the machining of the gear. In 
addition, two pits are observed on the surface. These pits could have resulted from 
several variables such as particles that were caught on the surface during meshing of 
the teeth, manufacturing flaws, surface wear, etc. The cuts and pits give an indication 
of the initial flaw size from which fatigue cracks may originate. 
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Figure B.2: Typical picture of flat, polished area on tooth #1 1 (400x). Photograph was 

taken near location B in Figure 6.4. 



Figure B.3: Machining grooves and pits on tooth surface (280x). 
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